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1.  Introduction 


Previously,  we  developed  a  method  to  allow  mammographic  differential  diagnosis  based 
upon  the  3-D  orientation  and  morphology  of  breast  calcifications.  This  method  used  a 
limited-view,  binary  reconstruction  technique.  In  clinical  trials,  it  was  shown  to  be  of 
value  in  instances  where  calcifications  are  associated  with  a  mass.  In  such  cases,  we 
could  distinguish  between  preferentially  peripherally  distributed  calcifications  that  are 
predominantly  benign  and  homogeneously  distributed  calcifications  that  are  more  likely 
to  be  malignant.  We  have  also  been  able  to  elucidate  the  linear  distribution  of 
calcifications  contained  within  a  ductal  system.  Unfortunately,  this  reconstruction  method 
does  not  allow  one  to  image  non-calcified  tissues,  or  relate  calcifications  to  the 
surrounding  tissue.  Thus,  in  the  present  study  we  have  developed  and  propose  to 
compare  alternative  methods  of  generating  3-D  images  of  the  breast,  namely  stereoscopy, 
linear  tomosynthesis,  limited  view  image  reconstruction  using  algebraic  reconstruction 
techniques,  and  micro-CT.  Such  methods  are  in  theory  capable  of  producing  3-D  images 
of  both  calcifications  and  the  surrounding  breast  tissue. 

It  is  our  purpose  to  develop  a  method  that  will  be  clinically  viable  in  terms  of  dose,  image 
quality  and  equipment  cost.  We  believe  that  these  proposed  developments  will  further 
enhance  the  3-D  imaging  and  evaluation  of  breast  cancer  by  allowing  the  radiologist  to 
view  calcifications  in  relation  to  the  surrounding  tissue,  and  to  allow  3-D  imaging  of  non- 
calcified  breast  tissues  at  doses  which  are  clinically  acceptable.  Stereoscopy  has  the  value 
of  providing  depth  perception  of  tissues  with  little  additional  dose,  however,  often  the 
small  angle  separating  the  views  is  insufficient  to  completely  determine  the  causes  of 
superposition.  Tomosynthesis  requires  more  views  and  potentially  higher  dose,  but 
provides  better  separation  of  tissues.  Artifacts  from  the  reconstruction  algorithms  can 
blur  synthetic  tomograms.  CT,  while  providing  the  best  3-D  images  requires  doses  that 
are  not  clinically  acceptable. 

To  date,  we  have  developed  an  imaging  system  capable  of  acquiring  each  of  the  images 
required  for  the  different  techniques.  We  are  in  the  process  of  quantifying  the 
performance  of  that  system,  and  preparing  to  acquire  specimen  images.  This  work  is 
reviewed  in  this  annual  report. 


5 


2.  Body 

2.1.  Summary  of  Work  Items 

It  is  useful  to  restate  the  work  items  listed  in  the  original  grant.  They  are  as  follows: 

Task  l:Develop  3-D  imaging  techniques  (Months  1-24) 

Subtask  la:  Stereoscopy  (Months  1-4) 

Subtask  lb:  Linear  Tomosynthesis  (Months  5-8) 

Subtask  lc:  Limited  View  Reconstruction  (Months  9-20) 

Subtask  Id:  Computed  Tomography  (Months  18-24) 

Task  2:Phantom  Development  (Months  1-24) 

Task  3:Evaluate  3-D  imaging  methods  with  phantoms  (Months  1-36) 

Task  4:Acquire  biopsy  specimen  image  datasets  (Months  13-36) 

Task  5:Observer  study  using  specimen  images  (Months  25-36) 

At  the  current  time,  tasks  1-3  are  ongoing,  while  tasks  4  and  5  have  not  yet  begun.  In  the 
following  report,  a  discussion  of  the  accomplishments  for  the  period  of  September  15,  1998  until 
March  2001  will  be  provided.  Due  to  extensions  and  delays  incurred  during  the  performance  of 
this  grant,  a  detailed  summary  of  the  timeline  of  events  affecting  this  grant  will  also  be  offered. 

2.2.  Development  and  Evaluation  of  3-D  Imaging  Techniques 

There  has  been  considerable  interest  in  imaging  the  breast  in  3-D.  This  has  included  numerous 
radiographic  methods  including  stereoscopy1,2,  tomosynthesis3,  limited-view  reconstruction  of 
calcifications4"6,  limited-view  tomography  reconstruction,  and  computed  tomography7.  There 
has  also  been  interest  in  3-D  ultrasound,  3-D  MRI,  and  other  potential  methods  of  imaging  the 
breast  in  3-D.  However,  to  date,  there  has  been  a  sizeable  gap  between  the  proposed  techniques 
and  clinical  feasibility.  To  achieve  clinical  feasibility,  it  is  necessary  to  consider  these 
techniques  on  the  basis  of  dose,  fundamental  imaging  physics,  and  technology.  This  grant  was 
designed  to  consider  4  such  methods  and  compare  them,  both  on  the  basis  of  physical 
performance  and  clinical  images. 

To  date,  we  have  had  a  number  of  difficulties  completing  this  grant.  These  difficulties  are 
described  in  detail  in  section  2.4.  However,  significant  changes  in  clinical  staff,  clinical 
demands  on  the  time  of  the  participants  of  this  grant,  and  the  sale  of  a  key  piece  of  equipment  by 
Thomas  Jefferson  University  Hospital  have  all  served  to  retard  progress  on  this  grant.  In  spite  of 
these  difficulties,  we  have  been  able  to  develop  alternative  solutions,  and  have  steadily  made 
progress  on  the  originally  proposed  work  items.  The  following  report  details  the  progress  to  date, 
discusses  the  difficulties  experienced,  and  discusses  solutions  to  these  problems. 

Image  Acquisition  Hardware 

The  grant  originally  proposed  to  use  two  pieces  of  hardware  for  the  completion  of  this  grant. 

The  first  was  a  prototype  device  built  in  Dr.  Maidment’s  laboratory.  It  consisted  of  an  Eikonix 
1412  linear  digital  camera,  coupled  to  an  x-ray  image  intensifier  (XRII)  via  lenses.  These  were 
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Figure  1  A  schematic  of  the  imaging  system.  The  system  consists  of  an  x-ray  tube,  a  system  of 
translation  and  rotation  stages  for  positioning  the  specimen,  an  x-ray  image  intensifier  (XRII),  lenses 
and  a  CCD  camera.  Pre-  and  post-object  collimation  is  not  shown. 

lenses.  These  were  mounted  on  an  optical  bench  with  a  Siemens  Bil50/30/50R  x-ray 
tube,  collimators,  and  a  specimen  holder  mounted  on  a  rotary  motion  stage.  A  Siemens 
Heliophos  5S  generator  supplied  the  x-ray  tube.  All  of  the  components  were  connected 
to  a  486-computer  running  Linux.  A  schematic  of  the  imaging  system  is  shown  in 
figure  1 .  This  system  was  intended  for  preliminary  investigations.  It  had  known  image 
quality  problems  that  would  prevent  it’s  use  through  the  entire  study.  Fortunately,  a 
second  system  existed  at  that  time  which  had  significantly  better  image  quality. 

The  second  system  consisted  of  a  Fischer  Mammotest  stereotactic  core  biopsy  system  and 
a  Fischer  MammoVision  digital  x-ray  detector  (Fischer  Imaging,  Denver,  CO).  We  have 
had  a  long-standing  relationship  with  the  Fischer  R&D  group,  and  Fischer  had  for  many 
years  provided  technical  support  for  this  work.  The  Fischer  system  was  supposed  to  be 
modified  to  include  a  computer-controlled  rotator  stage.  This  would  allow  us  to  mount  a 
specimen  holder  on  this  stage  and  use  the  Fischer  system  to  produce  tomographic  images. 
This  system  was  being  developed,  in  part,  because  of  its  physical  location  in  the  Breast 
Imaging  Center,  which  is  in  the  same  building  as  our  outpatient  breast  surgery  center. 

We  had  planned  to  interview  women  prior  to  needle  localization,  including  obtaining 
informed  consent.  Then,  following  their  surgery,  the  surgical  specimen  would  be  rapidly 
imaged  on  the  modified  Fischer  table.  Unfortunately,  the  hospital  sold  this  unit,  August 
10th,  1999,  just  8  weeks  after  I  had  finally  found  a  post-doctoral  researcher  to  fill  the 
position  necessary  to  fulfill  the  requirements  of  this  grant. 

Since  that  time,  we  have  made  a  concerted  effort  to  allow  this  grant  to  proceed.  We  have 
done  this  in  two  ways.  First,  we  have  altered  the  original  imaging  system  by  replacing 
the  image  intensifier  and  Eikonix  camera  with  two  other  detectors.  We  have  begun  to 
characterize  their  performance,  and  are  currently  awaiting  parts  to  assemble  the  final  unit 
and  begin  more  detailed  experiments  for  system  characterization,  optimization,  and  tissue 
imaging.  The  first  of  these  replacement  detectors  uses  a  DRC  (DRC  is  a  subsidiary  of 
Hologic,  Wilmington,  DE)  amorphous  selenium  active  matrix  detector.  The  second  uses 
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an  SMD  (SMD,  Fort  Collins,  CO)  1M30  CCD  camera  coupled  via  lenses  to  a  phosphor 
screen.  In  the  following  sections,  the  performance  of  the  original  Eikonix  detector  is 
compared  with  the  DRC  and  the  SMD  detectors.  Then,  a  summary  of  the  research  in  the 
four  imaging  modalities  is  presented 

Detector  Characterization 

In  each  of  the  designs  discussed  below,  a  Siemens  Heliophos  5S  x-ray  generator  and  a 
Siemens  Bil50/30/50R  x-ray  tube,  a  Parker  computer  controlled  rotary  stage,  and  a 
custom  specimen/phantom  holder  are  used.  These  form  the  basis  of  the  imaging  system 
being  built.  They  are  all  mounted  on  a  optical  breadboard  that  is  equipped  with 
appropriate  rails  and  mounting  hardware.  A  rotate  only  geometry  is  used  to  acquire 
images.  The  setup  allows  acquisition  of  either  tomographic  or  stereoscopic  images, 
depending  upon  the  acquisition  protocol  used,  by  varying  the  angle  of  acquisition  and  the 
dose  per  image.  Each  detector  is  capable  of  acquiring  2-D  images,  and  hence  3-D 
volume  reconstruction  is  also  possible. 

Eikonix  Detector 

The  first  detector  built  consisted  of  an  Eikonix  1412  digital  linear  camera  and  a  Siemens 
9”/6”  x-ray  image  intensifier  (XRII).  Because  this  was  a  linear  camera,  true  CT  images 
were  acquired  differently  than  the  other  methods.  In  the  limited  view  methods  (including 
stereoscopy),  2-D  images  were  acquired  by  scanning  the  detector  at  selected  angles. 
However,  due  to  time  constraints,  CT  images  were  made  in  two  different  ways.  Angular 
sub-sampling  allowed  us  to  acquire  up  to  200  images  for  3-D  reconstruction. 
Alternatively,  1-D  images  could  be  acquired  at  a  greater  number  of  angles,  but  only  a 
single  2-D  slice  of  the  object  would  be  reconstructed. 

This  camera  was  used  to  develop  the  image  reconstruction  algorithms  used  for  the  CT 
images  (discussed  in  detail  below).  In  the  original  grant  application,  the  CT  images  that 
were  shown,  were  of  very  high  contrast  objects,  and  were  acquired  as  1-D  samples  to 
produce  a  single  2-D  slice.  A  conventional  filtered  back-projection  reconstruction 
algorithm  was  used.  Since  that  time,  we  have  altered  our  image  acquisition  code  to  allow 
2-D  images  to  be  acquired,  and  now  use  a  modified  Feldkemp  algorithm  to  generate  a 
stack  of  2-D  slices,  which  are  rendered  as  3-D  volume  data  sets. 

The  Eikonix  camera  has  proven  to  be  difficult  to  work  with.  A  significant  problem 
occurred  that  related  to  the  IEEE-488  protocol  that  was  used  to  communicate  between  the 
camera  and  the  computer.  These  problems  took  approximately  2  months  to  address. 
When  complete,  the  image  quality  was  improved,  but  as  w'e  built  more  challenging 
phantoms,  we  found  another  problem.  This  related  to  the  use  of  an  XRII  in  the  detector 
system.  We  have  extensively  studied  this  problem.  The  problem  stems  from  subtle 
motion  of  the  image,  and  regional  varying  intensity  fluctuations  in  the  XRII.  The  result 
is  circular  and  cross  artifacts  illustrated  in  the  reconstructed  images  shown  in  figure  2. 

When  we  first  noticed  this  problem,  we  spoke  to  both  XRII  manufacturers  and  others  in 
the  microtomography  field.  The  manufacturers  disavowed  knowledge  of  the  problem, 
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Figure  2  Reconstructed  images  produced  with  the  Eikonix  detector.  The  left  image  is  a  wire  phantom 
designed  to  test  rotational  artifacts.  Note  that  the  wire  is  reconstructed  with  different  levels  of  success 
(they  should  look  like  dots,  not  crosses).  Also  note  the  circular  artifacts  due  to  the  instability  of  the  XRII. 
The  image  on  the  right  is  of  a  fresh  chicken  thighbone  (8.8  mm  x  7.4  mm),  clearly  showing  a  ring  of 
cortical  bone,  and  details  within  the  medullary  cavity.  This  image  still  demonstrates  subtle  circular 
artifacts  that  limit  the  detection  of  low  contrast  structures  in  the  bone. 

however  every  other  researcher  in  the  field  of  microtomography  who  uses  an  XRII  in 
his/her  imaging  system,  notices  similar  problems.  In  our  discussions,  it  became  evident 
that  most  people  solve  this  problem  by  using  such  systems  only  to  image  high  contrast 
objects.  Since,  we  had  planned  from  the  beginning  to  image  low  contrast  objects,  we 
needed  to  seek  a  different  solution.  Thus,  we  decided  to  change  detector  designs.  We 
came  to  this  decision  in  August,  2000. 


PRC  Detector 

At  the  time  of  our  decision  to  abandon  the  Eikonix  camera,  we  considered  two  possible 
designs.  The  first  was  based  upon  a  DRC  active  matrix  detector  that  we  had  in  our 
laboratory.  This  system  has  very  high  modulation  transfer  function  (MTF)  (see  figure  3), 
and  excellent  noise  power  spectra  (NPS)  and  detective  quantum  efficiency  (DQE).  These 
data  were  measured  in  our  laboratory.  In  support  of  this  work,  we  did  extensive 
modeling  and  experimentation.  This  work  was  reported  in  Medical  Physics  in  October 
2000,  a  reprint  of  which  is  included  with  this  report8.  The  system  has  been  used  to 
acquire  tomographic  images.  Examples  are  shown  in  figures  8  and  9  during  the 
discussion  of  computed  tomography 

The  great  strength  of  this  detector  is  the  image  quality.  We  intend  to  continue  to  use  it 
for  at  least  some  of  the  experiments  to  test  individual  3-D  acquisition  methods.  However, 
this  detector  is  very  slow,  producing  one  image  every  40  to  50  seconds.  Thus,  even  using 
angular-undersampling,  tomographic  images  produced  from  200  individual  views  take 
more  than  2  hours  to  acquire.  Thus,  we  had  to  look  for  a  different  detector  to  perform  our 
CT  research. 
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Ip/mm 

Figure  3  MTF  of  the  DRC  detector,  compared  with  the  ideal  MTF  of  a 
sine  function.  Note  that  a  small  amount  of  edge  enhancement  occurs  due 
to  the  use  of  an  amorphous  selenium  photoconductor 

SMD  Detector 

The  second  new  detector  design  was  based  upon  a  CCD  camera  that  was  coupled  with 
lenses  to  a  phosphor  screen.  This  design  allows  us  to  use  a  very  fast  camera.  The  camera 
we  chose  is  capable  of  imaging  at  30  frames  a  second.  This  allows  tomographic  images 
to  be  acquired  in  as  little  as  7  seconds  for  200  projections  to  34  seconds  for  1000 
projections.  However,  now  we  face  limitations  of  the  x-ray  generator,  which  will  extend 
that  time  to  a  few  minutes.  The  problem  with  this  design  is  that  it  is  only  ever  going  to 
be  useful  as  a  laboratory  system.  The  use  of  lenses  to  couple  light  from  a  phosphor 
screen  to  a  CCD  camera  is  quite  inefficient.  Our  calculations  have  shown  that  as  little  as 
2  light  quanta  per  x-ray  interaction  may  be  recorded.  We  have  shown  previously9’10  that 
this  inefficiency  can  result  in  needing  higher  doses  for  the  patients  (or  test  objects).  This 
is  not  considered  a  great  problem  for  the  experiments  currently  planned,  as  we  only 
intend  to  image  tissue  samples  and  inanimate  test  objects.  It  is  our  expectation  that  a 
clinical  system  would  use  one  of  the  newer  active  matrix  arrays  that  are  capable  of  10  to 
30  frames  per  second,  and  which  would  be  dose  efficient. 

We  obtained  a  demonstration  detector  from  the  manufacturer  during  the  period  of 
December  2000  to  February  2001.  During  that  time,  we  wrote  the  necessary  software  to 
control  the  camera  with  the  existing  x-ray  system.  We  performed  preliminary  testing  of 
the  camera,  and  were  satisfied  with  the  performance.  An  example  of  a  resolution  pattern 
is  shown  in  figure  4.  One  tomographic  image  was  also  obtained.  At  the  time  of  this 
report,  the  final  equipment  for  this  system  has  been  ordered.,  however,  it  will  not  be 
available  until  early  May  2001 .  In  the  interim,  we  continue  to  use  the  DRC  detector  to 
perform  the  Raleigh  discrimination  task  described  under  the  stereoscopy  section. 
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Figure  4  A  portion  of  a  resolution  phantom, 
showing  a  limiting  resolution  of  1 1.3  lp/mm, 
acquired  with  the  SMD  detector,  a  Lanex 
phosphor  screen,  and  85mm:55mm  relay  lenses 
giving  an  overall  field  of  view  of  20  mm.  A 
subsection  of  the  image  is  enhanced  for 
printing. 


Image  Technique  Optimization 

Stereoscopy 

Stereoscopy  is  the  process  by  which  two  views  of  a  scene,  obtained  at  slightly  different 
angles  (such  as  that  due  to  the  displacement  of  our  two  eyes)  provides  the  viewer  with  the 
perception  of  depth  (i.e.,  the  ability  to  discern  that  one  object  is  behind  another).  This 
process  can  be  applied  to  radiographic  imaging  to  achieve  a  similar  perception  of  depth. 
Originally,  we  proposed  to  perform  a  Raleigh  discrimination  test  to  determine  optimal 
angular  separation  for  radiographic  imaging.  Although  we  have  performed  a  extensive 
literature  search  to  help  design  the  experiment,  the  equipment  problems  described  above 
have  delayed  the  experiment.  The  experiment  is  currently  ongoing  using  the  DRC 
detector.  It  is  hoped  that  it  will  be  complete  prior  to  the  arrival  of  the  new  SMD  camera 
in  May  2001. 

Linear  Tomosynthesis 

We  have  also  performed  preliminary  experiments  with  regard  to  linear  tomosynthesis, 
using  the  Fischer  Mammovision  stereotactic  breast  imaging  system,  prior  to  its  sale.  For 
each  tomosynthetic  data  set,  fifteen  images  were  acquired  as  the  x-ray  tube  was  moved 
through  a  30  degrees  arc.  The  advantage  of  synthetic  tomography  over  conventional 
tomography  is  that  the  set  of  15  images  can  be  used  to  reconstruct  an  arbitrary  number  of 
planes,  while  in  conventional  tomography  each  image  plane  requires  an  additional  x-ray 
exposure 

.The  principles  of  tomosynthesis  are  illustrated  in  figure  5.  As  the  x-ray  focus  is  moved, 
the  imaging  plane  P  is  held  fixed.  To  perform  a  reconstruction,  each  x-ray  image  is 
viewed  as  a  gray-scale  function  gi  defined  on  a  region  Pr  of  the  imaging  plane  P,  where  i 
enumerates  the  x-ray  exposure.  For  each  plane  Q  parallel  to  P  in  which  reconstruction 
will  be  performed,  each  projection  position  of  the  x-ray  focus  defines  a  one-to-one 
correspondence  of  points  in  the  plane  Q  with  the  points  in  the  plane  P.  Explicitly,  for  the 
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i- th  position  of  the  x-ray  tube,  for  each  point  q  in  Q,  there  is  a  line  passing  though  the 
focus  and  q  which  meets  P  in  a  unique  point  p.  Thus  for  every  i  there  is  a  "projective" 
transformation  between  the  points  of  P  and  Q.  Further,  the  gray-scale  function  gi  on  Pr 
can  now  be  considered  as  a  function  g’;  defined  in  the  region  Q[  in  the  plane  Q,  and  the 
value  of  the  reconstructed  gray  scale  image  at  a  point  q  is  the  sum  of  all  the  g  ’j  which  are 
defined  at  that  point.  As  with  conventional  tomography,  for  objects  lying  within  the 
plane  Q,  the  functions  g )  add  coherently  to  produce  a  focused  image,  while  for  objects 
outside  the  plane  are  blurred. 

An  example  of  images  of  a  preliminary  tomosynthesis  phantom  is  shown  in  figures  6  and 
7.  The  phantom  consists  of  Lucite  spheres  and  cubes  contained  in  a  water-filled  Lucite 
box.  This  box  is  superimposed  upon  a  contrast-detail  phantom.  In  figure  6,  three 
projection  images  of  this  phantom  are  shown.  These  images  demonstrate  the  overlaying 
"clutter"  of  the  spheres  and  cubes  that  effectively  obscure  most  of  the  features  of  the 
contrast-detail  phantom.  Even  at  10  times  the  dose  this  image  was  acquired  at,  we  could 
not  visualize  more  than  5  objects  of  the  contrast-detail  phantom.  In  figure  7,  we  show 
three  reconstructions  of  the  phantom  at  different  heights  above  the  detector.  Two 
sections  contain  only  the  acrylic  spheres.  At  a  height  of  10  mm,  however,  the  contrast- 
detail  phantom  can  be  seen.  In  this  image,  24  contrast-detail  elements  are  visible.  The 
increase  in  detection  is  due  to  the  reduction  in  the  overlaying  clutter.  Thus,  reduction  of 
this  structural  noise  results  in  an  increase  in  the  effective  SNR  of  the  contrast-detail 
elements,  without  requiring  an  increase  in  dose. 


Figure  5  An  x-ray  tube  is  shown  at  three  locations  (1,2,3);  the  imaging  plane,  P,  is 
held  fixed.  An  image  of  the  letter  “A”,  located  in  a  plane  Q  is  projected  to  three 
different  locations  on  the  detector.  Note  that  the  individual  regions  Qi  overlap. 
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With  regard  to  our  process,  the  tomosynthesis  image  display  software  is  complete.  Image 
reconstruction  software  has  been  written,  but  further  refinements  are  necessary,  including 
removal  of  edge  artifacts  in  the  reconstructed  images.  Image  acquisition  optimization 
and  image  processing  optimization  also  must  still  be  completed. 


Figure  6  Three  source  images  of  a  tomosynthesis  phantom.  The  lead  BB’s  on  the  top  edge  act  as  fiducial  markers  for 
determining  the  angle  of  the  image.  The  phantom  consists  of  luctite  cubes  and  spheres  in  a  water  bath.  Attached  is  a 
contrast  detail  phantom.  The  low  contrast  elements  of  the  contrast  detail  phantom  are  not  readily  seen. 


Figure  7  The  tomosynthetic  images  of  the  phantom  shown  in  figure  XX.  The  phantom  is  shown  at  three  diffierent 
depths,  including  one  the  incorporates  the  contrast-detail  elements  (leftmost).  Numerous  elements  are  now  visible 


Limited-View  Reconstruction 

We  have  performed  an  extensive  literature  search  and  have  learnt  the  science  of 
tomographic  reconstruction,  however  little  experimental  work  has  occurred  with  regard 
to  this  work  item  to  date 

Computed  Tomography 

A  computed  microtomography  system  has  been  built  to  image  breast  specimens.  The 
system  is  currently  under  evaluation.  As  discussed  above,  the  system  has  been  built  so 
that  it  can  accommodate  two  different  types  of  detectors,  a  DRC  active  matrix  detector, 
and  one  based  on  a  phosphor  screen  optically-coupled  to  a  CCD  camera.  The  DRC 
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detector  is  complete,  while  the  CCD  based  detector  will  be  complete  in  May  2001  (a 
prototype  built  with  loaned  equipment  was  tested  between  December  2000  and  February 
2001). 

We  have  yet  to  finalize  the  data  acquisition  parameters.  We  are  evaluating  both 
conventional  and  fully  3-D  CT  image  acquisition  methods.  Experimentally,  we  have 
acquired  images  using  between  200  and  1000  projections,  rotating  over  an  angle  of  360°. 
Given  the  acquisition  geometry  used,  there  is  a  minimum  number  of  projections  that  are 
required  to  avoid  undersampling  the  projection  space.  This  number  depends  upon  many 
factors  including  the  object  size  and  detector  pitch.  However,  undersampling  generally 
only  causes  artifacts  in  very  highly  attenuating  objects  such  as  bone.  Thus,  given  the  low 
attenuation  of  breast  tissues,  undersampling  is  likely  to  be  less  of  a  problem.  This  is  one 
of  the  justifications  of  considering  limited  view  reconstruction  techniques  later  in  the 
grant. 

For  performing  simple  axial  reconstructions,  the  projection  data  are  reconstructed  using  a 
filtered  back-projection  algorithm  (RECLBL,  Donner  Laboratory).  A  Hanning  filter  is 
used  in  the  reconstruction.  Prior  to  reconstruction,  the  projection  data  are  corrected  for 
pixel-to-pixel  variations  in  the  detector  response,  fluctuations  in  x-ray  exposure  between 
projections,  and  error  in  the  center  of  rotation. 

Numerous  phantoms  have  been  constructed  and  imaged.  Two  are  shown  in  figure  8. 
Shown  are  a  tomographic  resolution  test  object  and  a  uniformity  test  object.  The 
resolution  test  object  has  5  rows  of  holes,  varying  in  size  from  1/16”  to  1/64”.  All  are 
clearly  visible.  The  uniformity  test  object  (right)  consists  of  a  Lucite  cylinder  1”  in 
diameter.  The  image  demonstrates  a  very  uniform  background  without  cupping  or  other 


Figure  8  Two  reconstructed  axial  images,  acquired  with  a  CT  scanner  based  on  the  DRC  detector.  The 
image  on  the  left  shows  a  resolution  phantom,  with  elements  of  size  1/16”  to  1/64”  all  visible.  On  the 
right  is  a  uniform  lucite  phantom,  clearly  showing  the  image  uniformity  and  low  noise  that  is  achievable. 
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artifacts.  These  images  show  the  excellent  potential  of  using  the  DRC  detector  to 
produce  tomographic  images. 


Shown  in  figure  9  are  some  preliminary  3-D  images  of  biological  tissue.  The  3-D  images 
are  rendered  in  two  different  ways,,  surface  rendering  and  volume  rendering.  Both  show 
exquisite  detail.  We  believe  that  these  are  the  first  CT  images  of  periplaneta  Americana. 
We  chose  to  image  this  species  due  to  it  ready  availability  and  due  to  the  size  of  the 
anatomic  structures.  Simply  stated,  this  represents  a  hard  imaging  task.  Attention 
should  be  paid  to  the  quality  of  the  reconstruction  of  the  leg  muscles.  Note  that  the 
muscle  fibers  inside  the  insect’s  exoskeleton  are  clearly  visible  and  discemable. 

Similarly,  part  of  the  exoskeleton  is  seen.  The  muscles  are  less  than  1  mm  in  diameter 
and  less  than  200  microns  separate  them  at  their  closest  point. 


Figure  9  Surface  rendered  (left),  and  volume  rendered  (right)  microtomography  images  of  species  periplaneta 
Americana  (common  American  cockroach),  acquired  with  the  DRC  detector.  The  insect  is  approximately  25 
mm  long.  The  muscle  fibers  and  the  exoskeleton  of  the  legs  are  clearly  visible. 


15 


2.3.  Discussion  and  Summary  of  Scientific  Results 

An  experimental  framework  for  performing  3-D  imaging  of  the  breast  has  been 
developed.  This  framework  consists  of  an  x-ray  generator,  x-ray  tube,  collimators,  and 
rotary  stage  for  holding  specimens  or  phantoms.  Each  of  the  components  is  under 
computer  control.  The  framework  is  designed  to  use  standard  optical  mounting  hardware 
and  an  optical  breadboard  to  allow  testing  of  many  different  detectors  and  acquisition 
protocols. 

We  have  developed  and  begun  to  characterize  two  detectors  that  will  replace  the 
detectors  originally  proposed  for  use  in  the  grant.  The  first  detector  is  an  active  matrix 
amorphous-selenium  device,  and  second  uses  a  phosphor  screen  optically  coupled  to  a 
CCD  camera.  The  first  device  has  excellent  dose  efficiency  (DQE)  and  resolution,  but  is 
too  slow  to  allow  easy  tomographic  imaging.  The  latter  device  again  has  excellent 
resolution,  but  lacks  the  dose  efficiency  of  the  flat  panel  detector.  This  detector  does, 
however,  have  the  temporal  response  needed  to  perform  CT. 

Initial  investigations  of  the  four  reconstruction  methods  have  been  delayed  due  to  the  sale 
of  the  original  detector  system  by  Thomas  Jefferson  University  Hospital.  However,  the 
new  detectors  will  both  be  working  soon,  and  work  on  the  reconstruction  methods  can 
then  proceed. 

Due  to  the  sale  of  the  Fischer  Mammotest  system,  it  will  be  necessary  to  alter  the  original 
work  plan.  The  new  detectors  are  not  physically  near  the  breast  imaging  and  breast 
surgery  centers.  As  such,  we  will  not  be  able  to  image  specimens  as  then  come  from 
surgery.  Instead,  we  will  image  mastectomy  and  cadaveric  breast  tissue. 
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2.4.  Discussion  of  Administrative  Issues 

At  this  time,  Dr.  Maidment  and  researchers  in  his  lab  are  effectively  working  towards 
completing  the  originally  proposed  work  items.  It  is  useful,  however,  to  discuss  the 
timeline  and  events  that  have  shaped  the  life  of  this  research  project.  When  first 
proposed,  the  work  was  to  be  performed  by  Andrew  Maidment  (PI),  Michael  Albert 
(research  assistant),  and  Emily  Conant  (clinical  collaborator).  Two  additional 
radiologists,  and  a  pathologist  were  included  to  provide  additional  clinical  assistance 
when  necessary.  Prior  to  beginning  the  grant  Dr.  Conant  left  Thomas  Jefferson 
University.  This  was  a  major  loss  as  Dr.  Conant,  among  all  of  the  radiologists  at  Thomas 
Jefferson  University  Hospital,  was  most  familiar  with  the  project  and  had  previously 
contributed  most  to  the  project.  We  have  never  since  had  a  radiologist  with  her  clarity  of 
insight  into  this  clinical  problem. 

Next,  since  the  departure  of  Dr.  Conant,  we  have  had  a  number  of  additional  people  leave 
the  TJUH  breast  center,  including  Dr.  Stephen  Feig,  Dr.  Dionne  Farria,  Dr.  Jane  Hughes, 
Dr.  Stephen  Lee,  Dr.  Steven  Nussbaum,  Dr.  Barbara  Cavanough,  and  others.  This  has 
had  a  negative  impact  on  all  of  the  research  being  performed  at  the  breast  center.  As 
previously  communicated  to  you,  Dr.  Catherine  Piccoli  assumed  Dr.  Conant’s  duties  after 
Dr.  Conant’s  departure,  and  Dr.  Farria  was  added  to  the  grant.  After  Dr.  Farria  left 
Jefferson,  her  position  on  the  grant  was  not  filled  by  another  person. 

Further,  due  to  clinical  duties,  Drs.  Maidment  and  Albert  were  expending  essentially 
100%  of  their  time  in  support  of  clinical  medical  physics  and  PACS  at  Jefferson.  The 
department  acknowledged  these  issues,  and  in  November  1998  hired  a  person  to  assist 
Dr.  Maidment.  No  salary  was  drawn  on  the  grant  until  November  1998.  In  spite  of  this 
assistance.  Dr.  Maidment  still  spends  at  least  40%  of  his  time  clinically.  In  addition,  Dr. 
Albert  had  sufficient  time  commitments  that  it  was  necessary  to  hire  a  post-doc  to  assist 
with  this  grant.  Finding  an  appropriate  post-doc  took  more  than  a  year.  However  in  June 
1999,  Dr.  Predrag  Bakic  began  working  on  this  project. 

Most  importantly,  however,  Thomas  Jefferson  University  Hospital  sold  a  key  piece  of 
research  equipment  on  August  10th,  1999  against  Dr.  Maidment’s  wishes.  This  device, 
the  Fischer  Mammotest/MammoVision  digital  x-ray  imaging  system  was  supposed  to 
have  been  used  in  virtually  all  of  the  experiments  planned  for  this  grant.  The  Radiology 
department  has  offered  to  help  by  providing  some  resources  that  would  offset  the 
capabilities  lost  with  the  sale  of  the  Fischer  system.  However,  as  of  this  date  (19  months 
later),  no  additional  resources  have  been  provided. 

As  a  result,  Dr.  Maidment  has  sought  other  methods  of  producing  tomographic  images. 
He  has  acquired  a  DRC  active  matrix  detector  and  a  related  computer  from  the 
manufacturer,  and  he  has  also  used  part  of  the  capital  equipment  budget  of  the  grant  to 
purchase  an  SMD  1M30  30-fps  12-bit  digital  CCD  camera,  and  Pentium-class  computer 
with  fast  image  acquisition  and  storage  capabilities.  During  a  trial  period  for  this  camera 
(December  2000  -  February  2001),  the  SMD  camera  was  extensively  tested  and  software 
to  operate  the  camera  was  written. 
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Thus,  at  this  time,  Dr.  Maidment  has  effectively  equipped  his  lab  with  the  equipment 
necessary  to  complete  the  work  in  this  grant.  Dr.  Maidment  is  currently  negotiating  with 
the  radiology  department  for  a  commitment  to  cover  the  salaries  of  participants  in  this 
grant  during  a  1  year  extension  to  this  grant,  which  we  will  submit  to  the  DOD  shortly. 

One  of  the  consequences  of  the  sale  of  the  Fischer  system  is  that  it  will  no  longer  be 
possible  to  acquire  images  of  tissue  specimens  of  women  undergoing  surgical  breast 
biopsies.  Pending  approval  of  the  program  office  of  the  DOD,  it  is  our  intention  to  use 
tissue  from  mastectomy  specimens  (after  pathological  examination),  or  to  use  cadaveric 
breast  tissue.  We  have  not  yet  sought  permission  from  either  the  IRB  or  the  other 
departments  involved  (surgery,  pathology,  and  anatomy).  These  negotiations  will  begin 
in  the  near  future,  pending  DOD  approval. 

In  spite  of  these  setbacks,  Dr.  Maidment  still  wishes  to  complete  this  research.  The 
scientific  questions  that  were  proposed  in  the  grant  proposal  are  still  relevant.  No  other 
research  group  has  tackled  these  questions.  The  following  work  remains:  optimization 
of  image  acquisition  techniques,  selection  of  appropriate  tissue  samples,  imaging  and 
evaluation  of  soft-tissue  images. 
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3.  Key  Research  Accomplishments 

The  following  is  a  list  of  key  research  accomplishments  resulting  from  this  work: 

•  Developed  a  generic  testbed  for  3-D  imaging  research,  consisting  of  x-ray  generator, 
rotary  stage,  specimen/phantom  holder,  and  detector  assembly. 

•  Developed  three  different  detectors 

o  the  original  was  based  on  an  Eikonix  linear  CCD  coupled  to  an  XR1I 
o  the  second  was  based  on  a  DRC  active  matrix  detector 
o  the  third  was  based  on  an  SMD  30  fps  CCD  camera 

•  Developed  stereotactic,  tomosynthetic,  and  volume  rendering  display  software 

•  Developed  theoretical  framework  for  measuring  and  calculating  the  physical  performance 
of  aliasing  and  non-aliasing  digital  x-ray  detectors. 

•  Begun  optimization  of  stereoscopy 

•  Begun  characterization  and  optimization  of  microtomography 
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4.  Reportable  Outcomes 

a)  Published  Manuscripts 

Michael  Albert,  and  Andrew  D.A.  Maidment.  Linear  Response  Theory  for 
Detectors  Consisting  of  Discrete  Arrays.  Medical  Physics,  27(10),  2417-2434, 
October  2000. 

b)  Abstracts  and  Presentations 

A. D.A.  Maidment,  3-D  Imaging  of  the  Breast.  6th  International  Cambridge 
Conference  on  Breast  Cancer  Screening.  Cambridge,  England.  April  14,  1999. 
(Invited  Presentation) 

A.D.A.  Maidment.  "3-D  Imaging  of  the  female  breast".  Imaging  2000, 
Stockholm,  Sweden,  June  29,  2000.  (Invited  Presentation) 

A.D.A.  Maidment,  P.  Bakic  and  M.  Albert.  "3-D  Digital  Mammography:  A 
Comparison  of  Image  Reconstruction  Methods".  DOD  Era  of  Hope,  Altanta,  GA, 
June  8-11,  2000. 

c)  Funding  Applications 

None 
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5.  Conclusions 


Initial  investigations  of  the  four  proposed  reconstruction  methods  have  been  delayed  due 
to  the  sale  of  the  original  detector  system  by  Thomas  Jefferson  University  Hospital.  To 
overcome  this  deficit,  an  experimental  framework  for  performing  3-D  imaging  of  the 
breast  has  been  developed.  This  framework  consists  of  an  x-ray  generator,  x-ray  tube, 
collimators,  and  rotary  stage  for  holding  specimens  or  phantoms.  Each  of  the 
components  is  under  computer  control.  The  framework  is  designed  to  use  standard 
optical  mounting  hardware  and  an  optical  breadboard  to  allow  testing  of  many  different 
detectors  and  acquisition  protocols. 

We  have  developed  and  begun  to  characterize  two  detectors  that  will  replace  the 
detectors  originally  proposed  for  use  in  the  grant.  The  first  detector  is  an  active  matrix 
amorphous-selenium  device,  and  second  uses  a  phosphor  screen  optically  coupled  to  a 
CCD  camera.  The  first  device  has  excellent  dose  efficiency  (DQE)  and  resolution,  but  is 
too  slow  to  allow  easy  tomographic  imaging.  The  latter  device  again  has  excellent 
resolution,  but  lacks  the  dose  efficiency  of  the  flat  panel  detector.  This  detector  does, 
however,  have  the  temporal  response  needed  to  perform  CT. 

Due  to  the  sale  of  the  Fischer  Mammotest  system,  it  will  be  necessary  to  alter  the  original 
work  plan.  The  new  detectors  are  not  physically  near  the  breast  imaging  and  breast 
surgery  centers.  As  such,  we  will  not  be  able  to  image  surgical  specimens.  Instead,  we 
will  image  mastectomy  and  cadaveric  breast  tissue. 

The  delays  introduced  by  the  sale  of  the  equipment,  and  personnel  issues  have 
sufficiently  delayed  the  completion  of  this  grant  that  a  no-cost  1  year  extension  will  be 
needed.  An  application  for  this  extension  will  be  submitted  shortly. 
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The  optical  transfer  function  (OTF)  and  the  noise  power  or  Wiener  spectrum  are  defined  for 
detectors  consisting  of  a  lattice  of  discrete  elements  with  the  assumptions  of  linear  response, 
Gaussian  statistics,  and  stationarity  under  the  discrete  group  of  translations  which  leave  the  lattice 
fixed.  For  the  idealized  classification  task  of  determining  the  presence  or  absence  of  a  signal  under 
signal  known  exactly/background  known  exactly  (SKE/BKE)  conditions,  the  Wiener  spectrum,  the 
OTF,  along  with  an  analog  of  the  gray-scale  transfer  characteristic,  determine  the  signal-to-noise 
ratio  (SNR),  which  quantifies  the  ability  of  an  ideal  observer  to  perform  this  task.  While  this  result 
is  similar  to  the  established  result  for  continuous  detectors,  such  as  screen-film  systems,  the  theory 
of  discrete  lattices  of  detectors  must  take  into  account  the  fact  that  the  lattice  only  supports  a 
bounded  but  (in  the  limit  of  a  detector  of  arbitrarily  great  extent)  continuous  range  of  frequencies. 
Incident  signals  with  higher  spatial  frequencies  appear  in  the  data  at  lower  aliased  frequencies,  and 
there  are  pairs  of  signals  which  are  not  distinguishable  by  the  detector  (the  SNR  vanishes  for  the 
task  of  distinguishing  such  signals).  Further,  the  SNR  will  in  general  change  if  the  signal  is  spatially 
displaced  by  a  fraction  of  the  lattice  spacing,  although  this  change  will  be  small  for  objects  larger 
than  a  single  pixel.  Some  of  the  trade-offs  involved  in  detectors  of  this  sort,  particularly  in  dealing 
with  signal  frequencies  above  those  supported  by  the  lattice,  are  studied  in  a  simple  model. 
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I.  INTRODUCTION 

The  importance  of  signal  detection  theory  in  quantifying  the 
performance  of  medical  imaging  systems  (x-ray  screen-film 
imaging  being  perhaps  the  best  example)  gives  impetus  to 
applying  the  same  techniques  to  the  digital  radiographic  im¬ 
aging  systems  which  are  now  coming  into  clinical  use.  As 
applied  to  screen-film  systems,  signal  detection  theory  re¬ 
quires  three  assumptions  to  be  at  least  approximately  ful¬ 
filled:  that  the  detector  responds  linearly  to  the  incoming 
signal,  and  is  both  stationary  and  homogeneous  (i.e.,  both  the 
detector  response  and  the  additive  noise  are  translationally 
invariant).  One  can  then  summarize  the  response  of  the  sys¬ 
tem  in  terms  of  the  gray-scale  transfer  characteristic,  the  op¬ 
tical  transfer  function  (OTF),  and  the  noise  power  or  Wiener 
spectrum. 

The  digital  x-ray  imaging  systems  which  are  now  appear¬ 
ing  generally  behave  as  a  lattice  of  discrete  detector  ele¬ 
ments.  Although  digital,  these  detectors  are  generally  oper¬ 
ated  under  conditions  such  that  the  effects  of  quantization  are 
negligible.  When  compared  to  screen-film  systems,  these  de¬ 
tectors  tend  to  be  linear  over  a  wider  range  of  exposures. 
Like  screen-film,  for  low-contrast  signals  the  noise  is  ap¬ 
proximately  additive  and  Gaussian.  However,  as  the  size  of 
the  imaging  elements  is  now  comparable  to  the  size  of  some 
of  the  smaller  objects  which  are  of  clinical  interest  (around 
0.1  mm),  these  detectors  are  not  strictly  homogeneous  in  that 
translations  by  a  fraction  of  the  lattice  spacing  result  in  the 
signal  being  recorded  in  a  different  manner.  As  these  devices 
generally  consist  of  a  regular  lattice  of  sensitive  elements, 
they  still  possess  a  symmetry  with  respect  to  a  discrete  group 


of  translations.  This  symmetry  is  approximate  due  to  the  fi¬ 
nite  extent  of  physical  detectors.  However,  as  in  the  theory 
of  screen-film  systems,  corrections  for  the  limited  extent  of 
the  detector  are  negligible  for  many  practical  applications. 
Thus,  one  can  apply  Fourier  techniques  to  put  the  signal 
detection  theory  of  such  devices  in  a  form  which  is  both 
tractable  and  similar  to  the  theory  of  screen-film  systems. 
Instead  of  using  a  continuous  Fourier  transform,  one  uses  a 
discrete  space  Fourier  transform,  which  recodes  the  data  ac¬ 
quired  by  the  detector  at  a  discrete  lattice  of  positions  in 
terms  of  a  bounded  and  continuous  range  of  spatial  frequen¬ 
cies. 

For  screen-film  systems,  the  OTF  diagonalizes  the  linear 
operator  which  relates  the  input  signal  to  the  output.  As  de¬ 
tailed  below,  for  discrete-array  detectors  the  effects  of  alias¬ 
ing  introduce  a  null  space,  different  for  each  device,  which 
prevents  this  operator  from  being  diagonalized  using  a  basis 
common  to  all  devices,  but  the  OTF  represents  the  operator 
in  a  basis  in  which  it  is  sparse  in  the  sense  that  all  terms 
vanish  except  those  between  input  and  output  spatial  fre¬ 
quencies  which  are  equal  or  aliased.  The  Wiener  spectrum  is 
the  discrete  space  Fourier  transform  of  the  discrete  autoco¬ 
variance  function,  and  thus  is  also  defined  in  the  region  of 
frequency  space  which  the  lattice  supports.  As  in  the  case  of 
continuous  detectors,  for  low-contrast  objects  (so  that  re¬ 
sponses  are  approximately  linear),  these  quantities  determine 
the  signal-to-noise  ratio  (SNR)  which  is  an  appropriate 
figure-of-merit  for  the  classification  task  of  discriminating 
between  the  presence  or  absence  of  an  exactly  known  signal 
against  an  exactly  known  background  (SKE/BKE). 
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Fig.  1.  (a)  a  rectangular  lattice,  (b)  The  reciprocal  lattice  of  (a),  (c)  A 
hexagonal  lattice,  (d)  the  reciprocal  lattice  of  (c).  Note  that  (a),  (b),  (c),  and 
(d)  represent  a  finite  region  of  a  lattice  which  covers  the  entire  plane,  (e)  A 
3X3  finite  rectangular  lattice,  (f)  The  circles  represent  the  frequencies  used 
in  the  finite  Fourier  transform  of  (e).  For  comparison,  a  unit  cell  of  the  full 
reciprocal  lattice  is  shown.  See  Sec.  Ill  for  details. 


unit  cell  (parallelogram  or  hexagon)  and  Fig.  1(d)  shows 
the  reciprocal  lattice  (note  that  v,  is  perpendicular  to  w2 
and  v2  is  perpendicular  to  w, ).  The  area  |A|  =  |viXv2|  of  a 
unit  cell  is  independent  of  the  choice  of  unit  cell,  since 
it  is  fixed  by  the  average  density  of  lattice  points  over  large 
regions.  The  area  of  the  unit  cell  of  the  reciprocal  lattice, 
M  =  |wi  Xw2|,  is  inversely  proportional  to  |A|,  as  can  be 
seen  by 


\A\\K\  = 


(vi)* 

(v2).c 


(Vl)y 

(V2)v 


det 


(wi), 

(Wl)y 


(w2).r\ 

(w 2)yl 


(19) 


/  V!  •  W[  V!  •  W-> 

det 

\  V2  *  Wj  v2-w2 


det 


1  0 
0  1 


=  1, 


(20) 


making  use  of  the  fact  that  the  determinant  of  a  product 
of  matrices  is  equal  to  the  product  of  the  determinants  and 
Eq.  (18). 

For  any  function  g(wt[, wt2)  of  the  lattice,  the  discrete 
space  Fourier  transform  is  defined1"'13  as 


|(f)  =  |A |  2  g(mi,m2)e  2lrir«i.»2‘f  (21) 

my  ,Wl2 

for  all  spatial  frequencies  f.  This  definition  is  equivalent  to 
evaluating  the  z-transform  on  the  unit  circle  in  the  complex 
plane.12'13  It  is  also  equivalent  to  the  Fourier  transform  of 
the  function  obtained  from  the  data  set  by  interpolation  with 


sine  functions  (e.g.,  Ref.  16,  p.  230).  Direct  calculation  from 
Eq.  (21)  gives 

g(f)  =  |(f+m,w1  +  m2w2),  (22) 

which  shows  that  g  is  periodic  in  Fourier  space  for  displace¬ 
ments  in  the  dual  lattice  and  one  need  only  consider  values 
of  g  on  one  unit  cell  of  this  lattice.  Any  frequency  f  outside 
of  this  unit  cell  is  an  alias  of  a  frequency  f  inside  the  cell, 
with  f-f  in  the  reciprocal  lattice.  Viewed  another  way,  the 
reciprocal  lattice  divides  points  in  the  frequency  plane  into 
equivalence  classes  of  points,  two  points  being  equivalent  if 
and  only  if  they  are  separated  by  a  vector  in  the  reciprocal 
lattice.  Any  unit  cell  will  contain  exactly  one  point  from  each 
equivalence  class  (except  for  boundaries),  and  knowledge  of 
g  on  the  unit  cell  determines  g  on  the  entire  plane.  Alterna¬ 
tively,  one  can  consider  g  as  being  defined  on  the  topological 
“quotient  space,”  a  torus,  just  as  one  can  consider  a  function 
on  the  real  line  with  period  2  tt  as  defined  on  the  unit  circle 
(Ref.  17,  p.  155). 

The  exponential  functions  in  the  discrete  Fourier  transfor¬ 
mation  satisfy  a  simple  orthogonality  condition 

1 1  d2fe-2witr^e2vit-r^.n2=\K\Smi .nSm^ni,  (23) 

where  K  is  the  region  corresponding  to  the  unit  cell  of  the 
reciprocal  lattice  in  the  frequency  plane  and  |AT|  is  the  area  of 
this  region,  thus  giving 

g(ml,m2)=  | fKd2f  i(D«2*if-,».-*  (24) 

as  the  inverse  transform.  The  complex  exponentials  form  a 
complete  set  of  orthogonal  functions,  so  that  any  appropriate 
periodic  function  of  frequency  f  can  be  represented  in  terms 
of  them.  The  completeness  can  also  be  expressed  in  terms  of 
a  comb  function  as 


2  e2™('-t'>  r’n^=\K\  2  (25) 

m  i  ^1 

where  the  equality  is  interpreted  in  terms  of  distributions  and 
the  sum  on  the  right-hand  side  is  over  the  frequencies  in  the 
reciprocal  lattice.8'16'17 

For  actual  finite  data  sets,  one  applies  the  finite  Fourier 
transformation.  The  discrete  space  Fourier  transformation 
can  be  interpreted  as  a  limit  of  the  finite  Fourier  transforma¬ 
tion  as  the  number  of  equally  spaced  points  in  the  data  set  is 
increased.  Specifically,  consider  a  bounded  subset  of  the 
{r„  ,,2}  such  as  M  of  the  form 

^g={r„|fJ Nx^nx<N\  ,N2^n2<N£,  (26) 

for  which  the  finite  Fourier  transform  and  its  inverse  are 
given  by 


gih,h)=  2  g(ni,n2)e 


-27n'f/ 


,n2 


(27) 


2  8Vi,h)e2 


1 1  ,l2 '  "  1  *n2. 


(28) 


where  \Ni=N'j-Ni.  The  points  in  the  Fourier  space  are 
given  by 
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h 

A  Ny 


Wi  + 


u 

A^W2 


and 


(29) 


&={titj2\Li*h<Ll>L2*h<L£, 


(30) 


where  the  U s  are  chosen  so  that  L\  —  L~N[  —Nt .  The  re¬ 
ciprocal  relationship  [Eq.  (18)]  between  the  basis  vectors 
{v,}  and  the  dual  basis  vectors  {w,}  gives 


*/|  ,/,'rn  |  ,n2 


/,  l2  \ 

a]v;w,+aa^W2)-(WiVi+"2V2) 


lxn{  l2n2 
AA^+ 


(31) 


which,  along  with  choosing  N>  =  -N  =N0I2  and  L,=  -L' 
=N0/2  for  N0  even,  produces  a  more  conventional  represen¬ 
tation  of  the  finite  Fourier  transform. 

As  the  number  of  data  points  ANXA.N2  increases,  the 
spacing  between  the  frequencies  f/  decreases,  so  that  in 
the  limit  the  data  points  on  the  lattice  extend  across  the  entire 
plane  and  the  frequency  values  fill  a  unit  cell  of  the  recipro¬ 
cal  lattice.  The  finite  sum  in  the  FFT  [Eq.  (27)]  approximates 
(with  a  factor  of  \A\)  the  infinite  sum  in  the  DFT  [Eq.  (21)], 
and  for  the  inverse  transform  the  sum  in  Eq.  (28)  (with  the 
introduction  of  a  factor  of  |A  |  [  AT|  =  1)  becomes 


<*> fft(^  i  >^2)  2 


1*1 


f AN  \AN2 


(|A|gFFr(^i  ^2)) 


(32) 


which  approximates  the  integral  used  in  the  inversion  of  the 
discrete  space  Fourier  transform,  Eq.  (24).  To  illustrate  this 
concept,  Fig.  1(e)  shows  a  small  rectangular  lattice  (corre¬ 
sponding  to  Ni=  —  1,  N\  —  2).  The  circles  in  Fig.  1(f)  repre¬ 
sent  the  corresponding  frequency  vectors  for  use  with  the 
finite  Fourier  transform.  The  box  shows  the  region  which 
would  correspond  to  a  unit  cell  of  the  reciprocal  lattice  if  the 
lattice  in  Fig.  1  (e)  were  extended  to  an  infinite  lattice.  If  the 
finite  lattice  shown  in  Fig.  1(e)  were  extended  (but  still  fi¬ 
nite),  the  corresponding  frequency  vectors  of  the  finite  Fou¬ 
rier  transform  would  fill  the  unit  cell  more  and  more  densely. 

It  should  be  noted  that  if  AN {  or  AN2  are  even,  then  some 
of  the  frequencies  at  which  the  finite  Fourier  transform  is 
defined  [shown  as  the  circles  in  Fig.  1(f)]  would  lie  on  the 
boundary  of  the  unit  cell,  and  such  frequencies  would  have 
aliases  which  also  lie  on  the  boundary.  For  example,  in  the 
square  lattice  considered  in  Fig.  1(f),  if  one  of  the  frequen¬ 
cies  at  which  the  finite  Fourier  transform  is  defined  fell  on 
the  edge  of  the  unit  cell,  an  alias  of  that  frequency  would  lie 
on  the  opposite  edge,  and  a  frequency  on  any  comer  would 
be  aliased  with  all  of  the  other  comers.  In  certain  sums  over 
frequency  components,  such  as  Eq.  (28),  it  is  useful  to  adopt 
the  convention  that  such  sums  include  exactly  one  represen¬ 
tative  from  each  class  of  aliased  frequencies,  so  that  frequen¬ 
cies  falling  on  the  boundary  of  the  unit  cell  are  not  counted 
multiple  times.  If  one  uses  the  “quotient  space”  point  of 


view,  this  follows  automatically  as  the  aliases  correspond  to 
a  single  point  in  the  quotient  space.  Alternatively,  one  might 
weigh  each  frequency  by  a  factor  (1/2  for  frequencies  lying 
on  edges  and  1/4  for  comers)  so  that  each  class  of  aliased 
frequencies  has  a  total  weight  of  1  (similar  to  counting  frac¬ 
tional  atoms  when  reckoning  the  number  of  atoms  in  a  unit 
cell  of  a  crystal). 

The  results  pertaining  to  Fourier  transformations  and  dual 
lattices  which  are  reviewed  in  this  section  have  direct  gener¬ 
alizations  to  any  number  of  dimensions,  but  as  the  statement 
of  the  results  for  arbitrary  finite  dimension  would  be  nota- 
tionally  cumbersome,  only  the  two-dimensional  results  have 
been  explicitly  stated.  For  notational  convenience,  let  m  rep¬ 
resent  the  ordered  pair  so  that  g(m)=g(m{,m2) 

^  rm=rml  ,m2>  and  similarly  for  k,  e.g.,  fk=fki  ^ 


IV.  TRANSFER  FUNCTION 

The  analog  of  the  optical  transfer  function,  which  relates 
the  response  of  the  detector  to  the  input  signal  in  frequency 
space,  can  now  be  defined.  The  input  signal  (/)  is  a  continu¬ 
ous  function  of  the  plane.  As  (/)  is  defined  relative  to  the 
“flat-field,”  it  is  reasonable  to  assume  that  (/)  has  compact 
support,  or  at  least  vanishes  sufficiently  quickly  at  infinity  to 
leave  the  quantities  considered  here  well  defined.  Thus  the 
Fourier  transform  (/)  is  a  continuous  function  of  the  entire 
frequency  plane.  The  data  D( rm)  are  well-defined  only  at  the 
discrete  lattice  points  rm,  so  that  the  discrete  space  Fourier 
transform  (D( f))  is  determined  by  its  values  in  one  unit  cell 
of  the  reciprocal  lattice.  Values  of  (D)  outside  of  the  first 
unit  cell  are  determined  by  the  periodicity  relative  to  the 
reciprocal  lattice  and  contain  no  new  information.  For  spatial 
frequencies  inside  the  first  unit  cell,  the  detector  responds  at 
the  same  frequency  as  the  input  signal.  For  frequencies  out¬ 
side  of  the  first  unit  cell,  the  detector  responds  at  an  aliased 
frequency,  so  it  is  impossible  to  uniquely  determine  the  input 
signal  without  additional  information,  although  it  will  be  ar¬ 
gued  in  later  sections  that  for  reasonable  tasks  this  is  not  a 
significant  problem. 

Each  point  on  the  detector  grid  is  assumed  to  respond 
linearly  to  the  incident  signal,  so  that  the  analog  of  Eq.  (1)  is 

(D(rm))  =  rJ  I  d2r'P(rm,r')(I(r')),  (33) 

where  P,  the  analog  of  the  point  spread  function,  represents 
the  response  of  the  detector  at  rm  to  x-ray  light  incident  at  r', 
and  T  is  a  constant  for  converting  x-ray  intensity  into  digital 
values,  generally  chosen  so  that  the  integral  of  P  with  respect 
to  r'  is  unity.  With  a  discrete  detector,  one  no  longer  has  full 
translational  invariance,  but  there  remains  an  invariance  un¬ 
der  translations  which  take  lattice  points  to  lattice  points, 
assuming  that  each  pixel  is  identical  except  for  position. 
Thus  we  can  write 

P(rm,r')  =  P(rm— r'),  (34) 

to  indicate  that  the  response  of  a  detector  element  to  an  input 
signal  depends  upon  the  displacement  of  the  detector  ele- 
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ment  from  the  region  to  which  the  signal  is  applied,  but  not 
upon  the  absolute  position  of  the  detector  element  or  the 
signal,  from  which  it  follows18,19  that 

(D(rm))  =  rJ  I  d2r'P( rm-r')</(r')>  (35) 

for  each  position  rm  of  the  sensitive  elements  on  the  lattice. 
While  the  data  D(rm)  are  only  available  at  the  lattice  points, 
the  convolution  can  be  calculated  at  any  point,  so  that 


r)  =  rj  |  d2r' P(r—r’)(I(rr)}, 

(36) 

J/(f)=rp(f)</(f)>, 

(37) 

serves  as  a  definition  of  ^(r)  for  any  position  r.  Although 
§)(r)  is  equal  to  the  data  (D(rJ)  at  the  lattice  points  where 
r=  rm ,  at  other  points  3r(r)  is  an  interpolation  which  will  not 
in  general  represent  a  physical  quantity,  although  it  is  some¬ 
times  useful  to  think  of  &{t)  as  the  response  of  a  virtual 
sensitive  element  added  to  the  detector  at  position  r  in  such 
a  manner  as  to  not  perturb  or  be  perturbed  by  the  other  ele¬ 
ments.  The  discrete  space  Fourier  transform  on  (Z)(rm))  can 
now  be  calculated  using  (D(rm))=^f(rm)  for  g(rm)  in  Eq. 
(21),  giving 

<D(I))  =  |A|2 

m 

=  |A|2  f  jKd2f,m')e2™m  (t'~r) 


<§(f+fk)=r2  r(f+fk)</(f+4)>, 

fk 


(38) 


which  follows  from  expressing  &  in  terms  of  its  Fourier 
transform  and  using  the  completeness  relationship  expressed 
in  Eq.  (25). 

Comparison  of  Eq.  (38)  with  its  screen-film  analog, 
Eq.  (2),  helps  to  clarify  the  interpretation  of  the  OTF,  T( f). 
The  spacings  in  the  discrete  lattice  introduce  new  length 
scales  which  occur  explicitly  in  the  summation  over  aliases. 
In  the  limit  of  a  very  finely  grained  lattice,  so  that  |A|— >0, 
the  spacing  of  the  reciprocal  lattice  points  gets  larger,  until 
only  the  one  unaliased  term  contributes  significantly  to  Eq. 
(38),  and  the  screen-film  case  is  recovered. 

When  frequencies  higher  than  those  supported  by  the  lat¬ 
tice  are  present  in  the  signal,  the  summation  in  Eq.  (38) 
introduces  “aliasing,”  that  is,  there  exist  multiple  spatial 
input  frequencies  whose  output  is  at  the  same  frequency  and 
are  thus  not  distinguishable.  For  example,  considering  a  one¬ 
dimensional  lattice  with  pixel-pitch  of  1  cm,  oscillations  at  a 
rate  of  0.5  cycles  per  cm  can  not  be  distinguished  from  os¬ 
cillations  at  a  rate  of  1.5  cycles  per  cm.  From  Eq.  (38),  two 
components  of  the  input  signal  generate  the  same  component 
of  the  output  signal  if  and  only  if  their  spatial  frequencies 
differ  by  an  element  tk  k  of  the  reciprocal  lattice. 

More  generally  for  any  lattice  there  are  frequencies  t0 
such  that  -  f0  is  an  alias  of  f0  (for  example,  if  f  and  -f  are 
on  opposite  boundaries  of  the  first  unit  cell  in  the  reciprocal 


lattice).  For  such  a  frequency  f0  [noting  that  T(f)  =  T*(  —  {) 
and  </( £))  =  </*(  — f»  for  real-valued  P( r)  and  </(r»]  it  is 
possible  to  choose  the  phase  of  (/(f0))  so  that 

r(fo)(i(fo))+T(-fo)<i(-fo))=0,  (39) 

showing  by  Eq.  (38)  that  a  sinusoidal  signal  concentrated 
at  frequency  f0  and  displaced  by  an  appropriate  offset  rela¬ 
tive  to  the  lattice  (as  determined  by  the  phase  of  (/(f„))) 
would  be  indistinguishable  from  the  flat-field  signal.  Return¬ 
ing  to  the  simple  one-dimensional  model  of  pixels  spaced  at 
1  cm,  this  result  means  that  for  some  displacement  relative 
to  the  lattice  the  input  of  a  sinusoidal  wave  of  frequency 
1  cycle/cm  would  give  vanishing  output.  If  the  detector  ele¬ 
ments  were  assumed  to  integrate  over  1  cm  intervals,  then 
the  output  vanishes  for  all  relative  phases  of  the  sinusoidal 
input  wave  and  the  lattice.  If,  alternatively,  the  detectors 
integrated  over  only  0.5  cm  regions  but  still  were  spaced  at 
1.0  cm  intervals,  then  the  sinusoidal  wave  would  have  van¬ 
ishing  output  only  when  the  nodes  of  the  sinusoid  fell  upon 
the  centers  of  the  0.5  cm  sensitive  regions  of  the  detectors 
and  would  otherwise  change  each  digital  value  by  a  phase- 
dependent  offset  from  the  flat-field  value. 

The  optical  transfer  function  has  been  written  in  terms  of 
a  Fourier  transform  using  complex  exponentials.  Since 
complex-valued  exponential  inputs  are  not  readily  available, 
it  is  necessary  to  ask  how  T  can  be  experimentally  measured. 
In  principle,  phantoms  machined  to  produce  sinusoidal  pat¬ 
terns  of  x-ray  intensity  could  be  used,  and  by  repeated  mea¬ 
surements  with  different  offsets  one  could  separate  the  posi¬ 
tive  and  negative  frequency  components.  A  more  practical 
method  is  the  well-known  slanted  edge  technique,20  21  in 
which  images  are  acquired  under  flat-field  conditions  except 
that  one  half-plane  of  the  detector  is  shielded  so  as  not  to 
receive  any  input  signal.  The  detector  response  D  as  a  func¬ 
tion  of  distance  from  the  edge  is  referred  to  as  the  edge- 
spread  function  ESF,  which  can  be  differentiated22  to  give 
the  line  spread  function,  LSF.  Alternatively,  by  providing  an 
appropriate  input  the  LSF  can  be  acquired  directly.23  The 
LSF  represents  integrals  through  the  PSF  along  lines  parallel 
to  the  edge,  so  that  by  acquiring  data  with  the  edge  at  mul¬ 
tiple  angles  one  obtains  the  radon  transform  of  the  PSF.  One 
can  reconstruct  the  PSF,  but  it  is  more  common  to  stop  after 
computing  the  Fourier  transform  of  the  ESF,  which  gives 
values  of  the  OTF  for  spatial  frequencies  f  which  are  normal 
to  the  edge.  For  discrete-array  detectors  it  is  desirable  that 
the  slope  of  the  edge  is  not  commensurate  with  the  lattice 
spacing  (for  example,  on  a  square  grid,  if  the  edge  is  not 
parallel  to  one  of  the  axes  and  does  not  have  a  slope  which  is 
a  ratio  of  small  whole  numbers  like  1/2  or  2/3).  When  this 
condition  is  satisfied,  for  a  given  region  of  interest  the  dis¬ 
tances  s  of  the  lattice  points  rm  from  the  edge  will  be  distrib¬ 
uted  sufficiently  densely  and  evenly  so  that  the  ESF  is  said  to 
be  “super-sampled,”  i.e.,  sampled  at  a  rate  significantly 
higher  than  the  reciprocal  of  the  lattice  spacing,  so  that  it  is 
possible  to  measure  values  of  the  OTF  for  input  frequencies 
beyond  those  supported  by  the  lattice. 

For  discrete-array  detectors,  rotational  symmetry  will 
generally  be  only  approximately  valid  at  low  spatial  frequen- 
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cies,  so  it  is  desirable  to  make  measurements  at  multiple 
angles  relative  to  the  lattice.  When  the  ESF^  at  a  given  angle 
6  is  acquired,  it  is  often  the  case  that  the  precise  position  of 
the  edge  relative  to  the  lattice  is  not  known,  so  that  one 
actually  acquires  data  for  ESF^z  +  z#),  where  z$  represents 
the  lack  of  knowledge  of  the  exact  position  of  the  edge. 
Upon  taking  the  Fourier  transform  of  the  ESF,  this  intro¬ 
duces  a  phase  uncertainty  of  the  form  e2m^Ze  into  the  value 
of  T(f).  While  this  phase  uncertainty  also  occurs  in  measure¬ 
ments  of  screen-film  systems,  for  discrete-array  systems 
summations  over  aliased  frequencies  generally  are  not  pos¬ 
sible  given  uncertainties  in  the  relative  phases  of  values  of  T 
at  different  spatial  frequencies.  In  general  one  can  remove 
this  phase  uncertainty  by  redefining  the  lattice  positions  to 
correspond  to  the  “centers-of-mass"  of  the  response  func¬ 
tions  of  the  sensitive  elements.  More  specifically,  if 

J  J  d2rP(r)=J  dz  ESF 0(z  +  zff)>O,  (40) 

then  it  is  possible  to  redefine  the  lattice  (by  a  shift)  so  that 
each  lattice  point  sits  at  the  center  of  mass  of  the  response 
function  of  the  associated  detector  element,  giving 

dxdy  xP(r)=J  J  dxdy  yP(r)  =  0.  (41) 

With  this  redefinition  of  the  lattice  position,  each  LSF  ac¬ 
quired  corresponds  to  a  radon  projection  of  the  PSF  (onto  a 
line  perpendicular  to  the  edge)  and  thus  the  center  of  mass  of 
each  LSF  should  be  at  the  origin.  This  corresponds  to  shift¬ 
ing  the  acquired  LSF  (adjusting  z&)  so  that 

J  dz  ESFfl(z)z  =  0  (42) 

for  each  angle. 

As  a  practical  matter  this  results  in  an  increase  in  the 
amount  of  data  it  is  desirable  to  report  for  a  given  detector.  If 
one  can  assume  an  inversion  symmetry,  i.e.,  P(r)  =  P(-r), 
then  the  imaginary  part  of  the  transfer  function  will  vanish 
identically,  so  that  only  the  real  part  need  be  reported.  The 
absolute  value  of  the  OTF,  traditionally  called  the  modu¬ 
lation  transfer  function  (MTF),  gives  enough  information 
to  calculate  quantities  such  as  the  spatial  average  of  SNR2 
(Sec.  VI),  but  does  not  give  enough  information  to  explore 
other  aspects  of  the  detector,  such  as  the  spatial  variation  of 
SNR2  as  the  test  object  is  moved  relative  to  the  lattice.  Re¬ 
searchers  should  also  note  that  with  the  slanted  edge  tech¬ 
nique,  when  combining  raster  lines  to  plot  the  edge  spread 
function,  the  independent  variable  of  interest  is  the  distance 
from  the  edge,  which  for  square  lattices  differs  from  the 
distance  along  a  raster  line  by  a  factor  of  the  cosine  of  the 
angle  between  the  raster  line  and  the  normal  to  the  edge. 
This  factor  becomes  significant  when  trying  to  measure  the 
transfer  function  at  angles  away  from  the  detector  axes. 
Based  upon  the  experience  of  the  authors,  one  can  generally 
measure  values  of  the  OTF  at  frequencies  several  times  the 
highest  frequency  supported  by  the  lattice.  One  is,  of  course, 
measuring  the  response  of  the  detector  at  low  frequency 
aliases  to  higher  frequency  input  signals.  Whether  the  pres- 
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ence  of  these  aliased  signals  in  the  output  is  desirable  will 
depend  upon  the  task  at  hand.  For  example,  it  might  be  de¬ 
sirable  to  detect  a  high-frequency  signal  even  if  one  can't 
distinguish  it  from  a  low-frequency  signal,  or  the  resulting 
ambiguity  might  be  unacceptable. 

The  question  of  what,  if  anything,  should  be  identified  as 
either  the  OTF  or  MTF  for  digital  systems  has  been  ad¬ 
dressed  in  several  ways  in  the  literature.  For  example, 
Dobbins24  discusses  the  “pre-sampled  OTF"  (OTFpre,  our 
T)  as  measured  via  the  LSF,23  but  then  emphasizes  the  fact 
that  the  response  to  an  input  signal  with  either  sinusoidal  or 
delta-function  spatial  variation  will  change  if  the  input  signal 
is  shifted  by  a  fraction  of  the  lattice  spacing.  This  depen¬ 
dence,  which  follows  from  Eq.  (38)  when  the  input  is  ex¬ 
panded  into  Fourier  components,  confounds  attempts  to  de¬ 
fine  the  MTF  either  in  terms  of  the  frequency  response  to  a 
single  delta  function  or  as  the  ratio  of  output-to-input  ampli¬ 
tude  for  a  sinusoid.  Dobbins  addresses  this  issue  by  defining 

OTFrf(f)  =  2  OTFpre(f+f|),  (43) 

fl 

and  defining  EMTF(f)  as  the  amplitude  of  the  detector  re¬ 
sponse  at  frequency  f  to  a  delta-function  input  averaged  over 
all  positions  of  the  delta  function.  Giger  and  Doi25  included 
such  a  summation  of  OTF  over  aliased  frequencies  in  their 
study  of  data  acquisition  and  display  for  digital  systems. 
Both  OTF^  and  EMTF  can  be  computed  in  terms  of  the  OTF, 
but  it  can  be  seen  that  neither  is  sufficient  for  calculating 
SNR.  Metz26  approaches  the  problem  in  essentially  the  same 
manner  as  discussed  in  this  paper,  and  indeed  Eqs.  (19)  and 
(30)  of  that  paper  essentially  give  our  Eq.  (38),  but  for  a 
slightly  more  specialized  case.  Metz  then  brings  up  the  point 
that  a  shift  by  a  fraction  of  the  lattice  spacing  in  the  input 
signal  does  not  result  in  a  simple  shift  in  the  output  data,  and 
concludes  that  “the  effect  is  accounted  for  mathematically, 
but  it  prevents  us  from  defining  a  unique  ‘transfer  function' 
of  the  sampling  process." 

Experimentally,  Sones  and  Barnes27  recognized  the  desir¬ 
ability  of  measuring  the  transfer  function  above  the  maxi¬ 
mum  frequency  supported  by  the  sampling  lattice  in  their 
work  with  a  digital  radiography  unit.  This  measurement  was 
performed  using  a  novel  technique  based  upon  a  phantom 
consisting  of  periodically  arranged  wires,  the  distance  be¬ 
tween  the  wires  chosen  to  be  incommensurate  with  the  dis¬ 
tance  between  samples  acquired  by  the  detector.  Fujita,  Doi, 
and  Giger28  measured  the  “pre-sampling  analog  MTF" 
above  the  maximum  frequency  supported  by  their  sampling 
lattice  via  a  slanted  slit  technique  and  recognized  that 
“knowledge  of  the  pre-sampling  analog  MTF  ...  will  be  use¬ 
ful  in  the  determination  of  signal-to-noise  ratio  (SNR)  [and] 
the  evaluation  of  digital  systems,"  a  statement  with  which 
we  heartily  agree. 

Working  from  a  complementary  theoretical  perspective, 
Barrett  et  al}  uses  the  “cross-talk"  matrix  to  address  the 
more  general  case  of  any  detector  whose  response  is  linear, 
then  proceeds  to  more  specialized  cases.  In  Barrett  et  al. ,  the 
input  to  the  system  is  defined  as  the  object  being  imaged 
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parameterized  in  terms  of  the  coefficients  of  its  three- 
dimensional  Fourier  series,  while  for  our  purposes  the  input 
is  the  x-ray  fluence  incident  on  the  detector.  For  projection 
radiography,  which  is  our  primary  interest,  the  incident  x-ray 
fluence  is  directly  related  to  the  integrated  attenuation  coef¬ 
ficient  of  the  object  along  rays  diverging  from  the  x-ray  fo¬ 
cus,  at  least  to  a  first  approximation.  As  our  goal  is  to  at¬ 
tempt  to  quantify  the  detector  response  independently  of 
other  technical  factors,  this  approximation  is  adequate.  Bar¬ 
rett  et  al.  is  concerned  with  detectors  which  may  have  rela¬ 
tively  few  sensitive  elements,  so  the  application  of  Fourier 
techniques  to  the  acquired  data  is  not  considered.  Barrett 
et  al.  applies  the  cross-talk  matrix  to  the  case  of  a  one¬ 
dimensional  array  of  detector  elements  with  aperture  size 
equal  to  the  element  spacing,  and  finds  that  the  cross-talk 
between  components  of  the  input  at  separate  frequencies  de¬ 
creases  as  the  length  of  the  array  is  increased,  so  long  as  the 
frequencies  are  not  aliases  of  each  other.  Thus  in  the  limit  of 
a  homogeneous  detector  of  infinite  extent  one  recovers  the 
fact  that  the  transfer  function  behaves  as  a  sparse  matrix,  in 
which  all  terms  vanish  except  those  on  the  diagonal  or  relat¬ 
ing  aliased  frequencies. 

In  order  to  use  Eq.  (38)  to  calculate  the  response  of  the 
detector  to  a  given  input,  it  would  be  necessary  to  know  the 
position  of  the  object  being  imaged  with  a  precision  finer 
than  the  lattice  spacing.  Strictly  speaking,  to  calculate  the 
response  in  either  the  discrete  or  continuous  case  requires 
that  the  input  be  “perfectly  known.”  However,  in  the  case  of 
a  continuous  detector,  a  shift  in  position  of  the  input  will 
result  in  a  corresponding  shift  in  position  of  the  output,  while 
for  a  discrete  detector  the  “shape”  of  the  output  would 
change.  In  many  cases,  such  as  predicting  the  detectability  of 
randomly  placed  objects,  one  would  need  to  calculate  for  an 
ensemble  of  objects  displaced  with  random  phases  relative  to 
the  lattice,  as  will  be  illustrated  below  in  calculating  the  SNR 
of  small  objects. 


V.  NOISE 

Individual  realizations  of  an  imaging  process  have  an  ir¬ 
reducible  variability  which  sets  a  fundamental  limit  on  how 
effectively  the  detector  can  distinguish  between  various  in¬ 
puts.  For  discrete-array  systems,  as  for  screen-film  systems, 
the  noise  can  be  quantified  in  terms  of  the  autocovariance 
function.  If  the  noise  is  additive  and  Gaussian,  then  the  au¬ 
tocovariance  matrix  completely  summarizes  the  stochastic 
process  which  generates  the  noise.  If  the  system  is  also  sta¬ 
tionary,  then  Fourier  techniques  can  be  used  to  define  the 
Wiener  spectrum.  i5 

The  discrete  autocovariance  function  is  given  by  '  ’ 

C(rm,rn)  =  <D(rm)D(rn)>,  (44) 

where  rm  and  rn  are  points  in  the  lattice  of  detectors,  the 
angled  brackets  represent  averaging  over  an  ensemble  of  flat- 
field  images,  and  as  discussed  above  (D(rm)}  =  0  in  the  ab¬ 
sence  of  a  signal.  Symmetry  under  interchange  of  positions 

C(rm,rn)  =  C(rn,rm)  (45) 
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is  an  immediate  result.  With  the  assumption  of  stationarity, 
the  autocovariance  depends  only  upon  the  displacement 
rm — rn,  so  we  can  write 

C(rm,r„)  =  C(rm-rn)  (46) 

without  ambiguity.  Note  that  the  difference  between  two 
vectors  corresponding  to  lattice  points  is  again  a  vector  cor¬ 
responding  to  a  lattice  point,  so  C  on  the  right-hand  side  of 
Eq.  (46)  is  defined  at  precisely  the  lattice  points. 

The  Wiener  spectrum  W( f)  is  defined  as  the  discrete 
space  Fourier  transform  [Eq.  (21)]  of  the  autocovariance 
function  C(rm).  As  with  any  discrete  space  Fourier  trans¬ 
form,  the  Wiener  spectrum  is  periodic  in  frequency  space 
[Eq.  (22)]  so  that  one  need  only  consider  the  values  of  W(  f) 
on  a  single  unit  cell  of  the  reciprocal  lattice.  It  is  noteworthy 
that  both  the  autocovariance  C  and  the  Wiener  spectrum  W 
are  real-valued  and  even.  As  with  screen-film  systems,  one 
considers  statistics  which  are  linear  functions  of  the  data,  so 
if  g(m{  ,m2)  is  a  set  of  real  (or  complex)  numbers  defined  on 
the  lattice  points,  one  defines 


0«.  =  |A|2  g(m)D(rJ. 
rm 

(47) 

The  variance  of  0g  (for  g  complex  valued,  the  sum 
variances  of  the  real  and  complex  parts)  is  given  by 

of  the 

Vm-(0g)  =  (0gd*) 

(48) 

=  |A|2((S  g(m)Z>(rj] ( 2  S*(n)D(r„)j 

) 

(49) 

=  |A|22  2  g(m)C(rn-rm)g*(n) 

m  n 

(50) 

in  terms  of  real  space.  Expressing  the  autocovariance  matrix 
as  the  inverse  discrete  space  Fourier  transform  [Eq.  (24)]  of 
the  Wiener  spectrum  one  obtains 

Var(^)  =  |A|22  2  g(m)  f  f  d2 f 

m  n  J  J  K 

X  W(  f)  e  2  m  f ' ( r"~  rm)g  *  ( n) 

(51) 

= f 

(52) 

where  the  second  step  follows  from  the  definition  of  the  dis¬ 
crete  space  Fourier  transform  [Eq.  (21)]. 

Thus  one  can  calculate  the  variance  of  a  statistic  6g , 
which  depends  in  a  linear  manner  upon  the  data,  using  either 
the  autocovariance  function  or  the  Wiener  spectrum.  Statis¬ 
tics  of  this  form,  for  g(m)  real-valued,  will  be  seen  to  cor¬ 
respond  to  decision  variables  of  ideal  observers  in  Sec.  VI. 
As  in  the  screen-film  case,  it  is  useful  to  consider  functions 
g(m)  corresponding  to  the  product  of  a  plane  wave  and  a 
windowing  function,  which  can  be  written  as 

gfo(m)  =  G(m)e2m'f‘’-r'",  (53) 

where  G(m)  is  a  real-valued  window  function  with  normal- 
ization 
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|A|2  G(m)G*(m)  =  1, 

m 


flf2fG(f)G*(f)=  1,  (54) 


where  the  two  normalizations  are  equivalent  by  Parseval’s 
theorem.  Applying  Eq.  (49)  and  Eq.  (52), 


Var( 0p)  =  |A|2(  2  G(m)e 


-2  it  if., 


mD(rJ 


(55) 


d2f\G(f-f0)\2W(f). 


(56) 


For  suitable  windowing  functions  G,  | G(f —  f0) |2  will  be 
strongly  peaked  near  f0  so  that  one  obtains  an  estimate  of  the 
Wiener  spectrum  at  the  specified  frequency,  W(f„).  In  par¬ 
ticular,  if  Grect(m)  is  chosen  as  l/(M lM1\A\)m  at  the  lattice 
points  rnGfO,...^!  —  1  ]X[0,...,M2-  1],  then 


M  \  —  \  M,- 1 


^rect  2  2 


mj=0  m2  =  0  \!lW  [M  2  A  \ 


D(vJe2*ifo"- 


|rrf-f  'i |2 _  sin2(A/17r(f— f0)~ Vi) 

11  °n  M,M2  sin2(7r(f— f0)- Vj) 

sin2(M27r(f-f0)-v2) 

X  sin2(ir(f— f„)  •  v2)  ’ 


(57) 


(58) 


which  explicitly  shows  that  for  this  choice  of  G, 
|G(f-f0)|2  is  strongly  peaked  near  (0 .  For  a  square  lattice 
with  conventional  choice  of  basis  vectors,  (f—f0)  •  Vj  =  (/* 
—  (A)*)  Ax,  where  fx — (f0)x  is  the  difference  in  the  x  com¬ 
ponents  of  the  frequencies  and  A*  is  the  lattice  spacing 
in  the  x  direction,  and  similarly  for  the  y  axis.  In  general, 
if  a  separable  window  is  chosen,  so  that  G(m) 
=  G1(wI)G2(m2),  then  G(f)  =  Gj(f* v1)G2(f- v2),  so  that 
one  can  make  use  of  the  variety  of  one-dimensional  windows 
which  have  been  studied.29 

Returning  to  the  case  of  a  general  lattice,  Eq.  (58)  shows 
that  for  this  particular  choice  of  window,  as  is  typical,  the 
estimate  of  W(f)  becomes  sharper  as  the  spatial  width  of  the 
window  increases,  so  that 


W(f)  =  C(f)=  lim  (W^^f)), 
M\M2^(X> 


|  A  | 


M  j-1  M2~  1 

2  2  D( rje-2*"-' 

m  |  =0  m t  =  0 


(59) 

l 

,  (60) 


where,  by  stationary,  any  M  { X  M2  region  of  the  detector 
lattice  will  serve.  Specializing  to  the  zero  frequency  case, 
f=0,  one  gets 


W(  0)=  lim  (M,M2|A|) 

.M2—*00 


X 


1 

M  {M 


Mx- 1  M2- I 


-2  2 

2  m i =0 


(61) 


which  is  the  discrete-array  version  of  Selwyn  granularity5 
(the  variance  in  the  average  digital  value  corresponds  to  the 


variance  in  the  spatially  averaged  optical  density  of  film). 
Comparing  Eq.  (61)  to  Eq.  (52),  one  can  interpret  Eq.  (61)  as 
the  statement  that  the  integrated  response  over  large  regions 
of  the  detector  depends  only  upon  the  low-frequency  com¬ 
ponents  of  the  Wiener  spectrum.  Viewed  spatially,  this  result 
means  that  the  digital  values  averaged  over  sufficiently  large 
disjoint  regions  are  approximately  independent,  so  that  the 
variance  of  the  average  over  N  large  subregions  scales  with 
l/Vocl/M,M2. 

As  with  the  OTF,  the  results  of  the  screen-film  theory 
appear  as  a  limiting  case  for  sufficiently  fine  lattices.  Writing 
Eqs.  (59)  and  (60)  as 


1 

^m{m2\a\ 

1  Mo-  1 

*( 

2  2 

wi|=  0  m 2  =  0 

(62) 


the  summations  become  approximations  of  the  integrals  in 
Eq.  (6). 

The  discrete  autocovariance  [Eq.  (44)],  the  definition  of 
the  Wiener  spectrum  as  the  discrete  space  Fourier  transform 
of  the  autocovariance,  and  the  use  of  Fourier  components  of 
flat-field  images  to  estimate  the  NPS  [Eq.  (59)]  have  oc¬ 
curred  in  several  places  in  the  medical  physics  literature,13,30 
but  historically  these  results  seem  to  have  been  considered 
less  than  satisfactory  from  a  theoretical  point  of  view.  For 
example,  Cunningham31  stated  that  while  “[i]t  is  tempting  to 
write  out  the  NPS  of  [the  sampled  digital  signal],  but  strictly 
speaking  this  violates  the  shift-invariance  assumption  since 
[the  data]  is  sampled  and  is  therefore  not  shift  invariant/ ' 
More  recently,  Cunningham,32  in  analyzing  the  concept  of 
NPS  in  terms  of  cyclostationary15'33  random  processes,  de¬ 
fines  WM^M^  [Eq.  (60)]  as  “a  working  definition  of  the  digi¬ 
tal  NPS/'  As  detailed  in  Sec.  VI,  the  NPS,  as  defined  here,  is 
precisely  the  noise  which  sets  the  detection-theoretic  limits 
on  the  use  of  the  detector.  In  the  detection-theoretic  approach 
of  Barrett  et  al.}  the  Fisher  information  matrix  relates  the 
detector  noise  back  into  uncertainties  in  the  estimates  of  the 
Fourier  coefficients  of  the  object  being  imaged.  This  has  the 
advantage  that  it  removes  the  fundamentally  arbitrary  choice 
of  scale  in  using  digital  values,  but  if  aliased  frequencies 
become  important  the  Fisher  information  matrix  becomes 
singular  so  that  the  inversion  of  this  matrix  is  problematic. 

The  definition  of  NPS  given  here  is  intended  to  be  opera¬ 
tional  in  the  sense  that  it  is  defined  in  a  manner  which  can 
be  implemented  using  the  experimentally  available  digital 
values.  For  the  purposes  of  understanding  the  sources  of 
noise  in  detectors,  it  may  be  useful  to  consider  the  noise 
in  the  “presampled’'  signal,  and  for  some  detectors  this  pre¬ 
sampled  signal  might  be  experimentally  accessible.  For  ex¬ 
ample,  in  a  detector  based  on  a  phosphor  screen  coupled  with 
a  lens  to  a  charge-coupled  device  (CCD)  camera,  one  could 
do  experiments  in  which  the  camera  is  replaced  by  a  photo¬ 
graphic  film.  For  some  devices,  such  as  TFT  arrays  using 
direct  conversion  mechanisms,  the  meaning  of  the  presa- 
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mpled  signal  is  less  clear  as  removal  or  refinement  of  the 
sampling  array  is  likely  to  change  the  electric  fields  respon¬ 
sible  for  charge  collection. 

As  reviewed  by  Wagner  and  Sandrik,30  the  calculation  of 
the  NPS  can  be  implemented  in  several  ways.  One  method  is 
to  estimate  the  autocovariance  function  [Eqs.  (44)  and  (46)] 
using  pairs  of  points  in  one  or  (preferably)  more  images,  and 
then  performing  the  Fourier  transform  to  give  the  NPS.  Al¬ 
ternatively,  the  variance  in  the  Fourier  components  is  used, 
as  in  Eq.  (55).  If  G  is  chosen  as  a  rectangular  window,  then 
Eq.  (55)  reduces  to  Eq.  (59),  so  that  ( [Eq.  (60)]  is 
used  as  an  estimate  of  W(f).  In  principle  the  frequency  f  is  a 
continuous  variable,  but  the  spread  of  |G(f)|2  limits  the  reso¬ 
lution  in  frequency  space  [by  Eq.  (56)]  and  this  spread  is 
inversely  proportional  to  the  size  of  the  spatial  region  and  on 
the  order  of  \K\IMyM2 ■  Given  this  resolution,  it  is  reason¬ 
able  to  calculate  the  NPS  at  MXM2  frequencies  spaced 
evenly  in  the  unit  cell  K  in  frequency  space.  Thus,  the  tech¬ 
niques  commonly  in  use  by  experimenters  give  precisely  the 
quantities  of  interest  from  our  current  theoretical  point  of 
view,  although  the  use  of  windows  other  than  the  rectangular 
window  might  be  of  interest  to  obtain  better  frequency  reso¬ 
lution. 

Generally,  frequency  resolution  is  not  a  limiting  factor  in 
estimating  the  Wiener  spectrum,  and  the  NPS  estimated  by 
(WM  Mft))  is  subjected  to  further  smoothing.  From  Eq.  (57) 
it  is  seen  that  <VFWim2( f„»  is  the  variance  in  the  random 
variable  0rect ,  and  as  the  region  of  interest  used  in  the  cal¬ 
culation  is  made  larger,  the  variance  in  0rect  tends  to  W(f) 
which  will  be  nonzero  in  general.  Because  the  variance  of 
6kci  does  not  vanish,  neither  will  the  variance  in  |0rect|2,  so 
the  variance  in  does  not  converge  to  zero  as 

M  | ,  M  2  ^  "c  ■  As  the  region  of  interest  is  made  larger,  one 
gains  in  spectral  resolution  but  not  precision,  and  this  repre¬ 
sents  an  unavoidable  trade-off.30'34  One  can  only  decrease 
the  uncertainty  in  the  estimates  of  the  Wiener  spectra  by 
averaging  estimates  of  14(f)  from  several  different  regions  of 
interest.  Of  course,  for  the  purposes  of  analysis  one  could 
divide  a  large  region  into  several  smaller  regions,  and  the 
averaged  value  of  estimates  of  14(f)  would  then  have  less 
uncertainty,  but  the  spectral  blur  would  be  increased.  Since  it 
is  often  inconvenient  to  obtain  sufficiently  many  flat-field 
images  to  make  the  standard  error  in  the  estimates  of  W(f)  at 
individual  frequencies  small,  researchers  often  opt  for 
smoothing  the  experimental  spectrum. 

VI.  KNOWN  SIGNAL  DETECTION 

Having  addressed  the  issues  of  OTF  and  Wiener  spec¬ 
trum,  it  is  now  possible  to  use  the  signal-to-noise  ratio  (SNR) 
to  quantify  the  ability  of  the  detector  to  perform  SKE/BKE 
tasks.  First,  however,  it  is  useful  to  briefly  review  the  mean¬ 
ing  of  the  SNR  in  terms  of  an  ideal9’35  observer  working  with 
Gaussian  statistics.  The  ideal  observer  is  challenged  with  de¬ 
ciding  between  two  hypotheses  based  upon  a  given  set  of 
data.  In  the  current  context,  these  data  consist  of  the  digital 
values  obtained  from  the  detector,  and  for  the  moment 
we  will  restrict  the  observer  to  knowledge  of  only  a  finite 


region  of  the  detector,  corresponding  to  indexes  me. 

This  observer  works 
under  the  assumption  that  given  hypothesis  Hl,  correspond¬ 
ing  to  an  expected  input  signal  (/( r»,  and  an  expected  data 
set  (£>(rm))j,  the  probability  density  function  describing  the 
expected  range  and  frequency  of  observed  data  sets  is  Gauss¬ 
ian.  This  Gaussian  distribution  in  (M[- M i)X(M[- M2) 
=  AM[AM2  dimensions,  one  dimension  for  each  detector 
element  available  to  the  observer,  can  be  written  explicitly, 
but  to  make  the  formulas  somewhat  less  cumbersome  we 
use  the  following  notation:  Xm=D(rm),  (Xm)i=(D(rm))I, 
(Xjn=(D(rJ)u,  and  {Xm}={D(  rm)|me^}  is  a 
A  M  j  A  M 2  -dimensional  vector  in  the  space  of  all  possible 
data  values  for  the  detector  elements  in  region  The  prob¬ 
ability  distribution  which  governs  the  frequency  with  which 
particular  data  sets  will  be  obtained  under  hypothesis  Hx  is 
given  by 


Pi({Xn})=N0e  l^™n*  /^Xm  1^mn^n 


where  the  normalization  factor  is  given  by 


(63) 


Nn 


J  \(AM1AM2)/2  j 

2  tt)  \/det  C 


(64) 


The  matrix  Cmn  is  the  autocovariance  function  C(rm,rn)  of 
Sec.  V  restricted  to  the  range  m,neJ.  The  fact  that  m  and 
n  are  double  indices,  e.g.,  m  stands  for  is  not  a 

problem  from  the  theoretical  point  of  view,  and  in  principle 
for  a  numerical  calculation  one  could  simply  choose  a  con¬ 
venient  one-to-one  pairing  of  the  double  indices  mx,m2 
e.yM  with  the  integers  l,...,AMiAM2  so  that  C  would  be 
indexed  in  a  more  customary  manner.  Under  hypothesis  Hu , 
the  range  and  frequency  of  observed  data  sets  will  be  gov¬ 
erned  by  a  Gaussian  probability  density  Pn,  this  time  con¬ 
centrated  around  (X)n.  The  restricted  covariance  matrix,  C, 
occurring  in  both  cases,  will  be  the  same  under  the  assump¬ 
tion  that  the  noise  is  additive. 

Returning  to  the  question  of  how  to  decide  between  hy¬ 
pothesis  Hx  and  hypothesis  Hu ,  if  for  a  given  instance  of  the 
experiment  a  data  set  {Xm}  =  {P>(  rm)|me^^}  is  obtained 
such  that  Pn({X J)  is  relatively  large  and  /^({XjO  is  rela¬ 
tively  small,  it  would  generally  be  reasonable  to  favor  Hu. 
Thus  the  ideal  observer’s  decision  rule  based  on  the  likeli¬ 
hood  ratio  Pn/Pj,  as  discussed  below,  is  intuitively  reason¬ 
able. 

The  ideal  observer  attempts  to  minimize  the  expected 
cost9,36  given  knowledge  of  the  cost  of  misclassification  un¬ 
der  either  hypothesis  and  the  a  priori  probabilities  associated 
with  each  hypothesis, 

(Cost)  =  P  ( Hi)  P  ( Chll|l)  C  i^n+  P(/7n)^(Chl|n)Cn^i , 


where  in  the  first  term  P(HX)  is  the  a  priori  probability  of  the 
state  corresponding  to  hypothesis  Hy  being  true,  P(ChIl|l)  is 
the  probability  of  mistakenly  choosing  hypothesis  Ha  when 
hypothesis  Hj  is  correct,  C j_jj  is  the  cost  associated  with  this 
error,  and  similarly  for  the  second  term.  Given  a  region 
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of  the  AM  { AM  2  dimensional  data  space  and  the  decision 
rule  that,  if  the  observed  data  {D(rm)|me.^}  are  in  7?n 
then  the  observer  rules  in  favor  of  hypothesis  Hu  and  other¬ 
wise  in  favor  of  Hl,  then  the  probability  of  mistakenly  fa¬ 
voring  hypothesis  Hn  when  Hi  is  correct  is 

P( ChIl|I)=  f  f  ...  f  dAM^-{xm}Pl({xm}),  (66) 
J  J  Jru 

and,  as  under  either  hypothesis  the  total  probability  must  be 
unity, 

P(Chl|H)  =  1  -  f  f  ...fR  dAM‘AM'-{XjPn({Xj)  (67) 

gives  the  probability  of  making  the  error  in  the  other  direc¬ 
tion.  Combining  Eqs.  (65)-(67), 

(Cost  }  =  P(Hn)Cn^ 

+  JI -\r  rfAM>AM2{Xm}DC({Xm}),  (68) 

where 

DC({Xm})  =  P(//I)C1^nPI({Xnl}) 


2  5./?(rm)C(rm,r„)  =  ((0(rn))n-(D(rn))I),  (74) 

me../# 

for  all  ne,/^.  On  physical  grounds,  the  values  of  the  mask 
function  g  ^(rm)  will  be  significant  only  in  the  region  near 
where  (7(r))n— (/(r))j  is  nonzero.  Further  from  this  region, 
the  values  of  g^(  rm)  will  tend  to  zero,  so  that  for  suffi¬ 
ciently  large  AMXAM2  the  ability  of  the  detector  to  dis¬ 
criminate  between  the  two  hypotheses  should  not  depend 
upon  the  exact  value  of  AM  x  AM  2 .  In  that  limit,  the  efficacy 
of  the  detector  for  the  SKE/BKE  task  should  be  set  by  the 
linear  statistic  6g  for  the  ideal  observer’s  mask  function  gj. 
This  mask  function  is  defined  implicitly  by 


\M  2  5i(rm)C(rni,rn)  =  ((Z)(rn))n-(D(rn))I),  (75) 

me../# 


where  a  factor  of  |  A  |  is  introduced  to  simplify  the  form  of 
the  solution  which  in  the  Fourier  domain  is  given  by 


U  f)= 


w(  f) 


(76) 


-P(Hn)Cn_iPn({Xm})  (69) 


is  the  differential  cost  which,  if  the  experiment  were  repeated 
sufficiently  often,  would  be  attributed  to  those  experiments 
which  gave  data  {D{ rm)|me  .#?}.  Clearly  the  expected  cost 
given  by  Eq.  (68)  is  minimized  by  choosing  the  region  Rn  to 
be  precisely  the  region  where  the  differential  cost  DC  is 
negative,  so  that  the  ideal  observer’s  decision  rule  is  to 
choose  hypothesis  Hu  if  and  only  if  the  likelihood  ratio 


Pn({Q(rm)jme./^}) 

P\({D(rJ\m^./6}) 

exceeds  the  threshold  value 

P(//i)Ct_n 
0  P(Hn)C^{ 


(70) 


(71) 


By  adjusting  the  operating  point  A0  one  makes  trade-offs  in 
the  rates  of  the  two  possible  error  types,  as  can  be  shown 
graphically  in  terms  of  receiver  operator  curves  (ROC 
analysis).37  Equivalently  one  can  place  the  cutoff  on  log  A, 
and  from  Eq.  (63), 


log  A=  2  ((Z)(rJ)n-(D(rm)>I)(C-,)mn7)(rn) 

m,ns.  /# 

+  const,  (72) 

where  the  constant  term  does  not  depend  upon  the  observed 
data.  Thus  an  ideal  observer,  viewing  a  finite  region  .  /#  of 
the  detector  array,  uses  a  linear  statistic  &  #  defined  by 

0//=  2  g./Arm)D(rm),  (73) 

me./# 

where  g  is  given  implicitly  by 


\({Ht+  fk))n-  </(f+  fk))i)  T(t+  fk) 

r— - 

w(  f) 


(77) 


The  statistic  6g  is  itself  a  Gaussian  variable  whose  variance 
can  be  computed  using  Eq.  (52),  so  that 


SNRf=r2 


sfk(A/(f+fk))r(f+fk)|2 

W) 


(78) 


gives  the  SNR  corresponding  to  the  use  of  the  statistic,  as 
defined  in  Eq.  (8).  Thus  the  limiting  case  of  a  detector  array 
of  infinite  extent  is  well  defined,  for  pixels  “far  away”  from 
the  region  of  interest  do  not  significantly  contribute  to  the 
decision.  Physically,  it  is  clear  that  the  “tails”  of  the  PSF 
and  autocovariance  functions  set  the  relevant  scale  by  which 
distance  from  the  edge  of  the  array  is  measured,  so  that  when 
the  projected  images  of  objects  appear  at  a  distance  from  the 
boundary  of  several  times  the  lengths  of  these  tails  the  de¬ 
tector  can  be  treated  as  essentially  infinite  and  Eq.  (78)  is 
valid. 

It  is  acknowledged  that  there  are  mathematical  subtleties 
related  to  a  truly  infinite  detector  which  are  not  addressed 
here.  For  example,38  the  data  set  for  such  a  detector  would 
represent  an  infinite  set  of  random  variables,  so  it  is  not 
possible  to  write  down  a  probability  density  distribution  like 
Eq.  (63)  in  the  infinite  case.  The  nature  of  the  physical  limit 
is  sufficiently  clear  that  a  study  of  these  mathematical  subtle¬ 
ties  could  not  change  the  results.  In  any  case,  the  fact  that  the 
linear  statistic  0g  with  g  —  gj  gives  the  optimal  SNR  of  any 
linear  statistic  can  be  proven  directly.  More  precisely,  if 
g(rm)  is  used  to  define  a  linear  statistic  6g ,  then  letting 
A  D  =  (D)n— (Z>)„ 
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1 1  j^2fg(f)AD(f)  2=j  f  J^d2f(g(f)JW)) 

J^)2 

W W)1 
JVf|s(f)l2w(f) 
r  r  I A  ml 

J  J  K  W(()  ’ 


(79) 


where  the  second  step  is  an  application  of  the  Schwarz  in¬ 
equality.  Dividing  both  sides  of  Eq.  (79)  by  the  first  factor  on 
the  right,  one  obtains 


\uKd2iSmm\2^  f  f  „>p(f)i2 

JJV2f|g(f)|2mf)  J  Jk  W(t) 


(80) 


where  the  quantity  on  the  left  is  the  SNR2  for  the  statistic  6g 
[Eqs.  (8)  and  (52)]  and  the  quantity  on  the  right,  proven  to  be 
larger,  is  the  SNR2  of  the  ideal  observer  as  given  by  Eq.  (78) 
[with  Eq.  (38)]. 

As  a  slightly  less  subtle  point,  the  construction  of  the 
ideal  observer  involves  dividing  by  W(f),  which  is  problem¬ 
atic  if  W(f)  =  0  at  some  frequency.  For  physical  detectors, 
the  Wiener  spectrum  never  vanishes  as  there  is  always  some 
residual  noise.  Even  for  highly  idealized  detectors,  the 
Wiener  spectrum  must  reflect  the  noise  in  the  incident  x-ray 
fluence  so  that  it  can  only  disappear  at  frequencies  where  the 
OTF  vanishes,  and  at  these  frequencies  the  Wiener  spectrum 
will  vanish  no  faster  than  OTF2(f)  (discussed  in  more  detail 
in  the  next  section),  so  that  even  in  this  case  the  SNR  as 
given  by  Eq.  (78)  is  a  well-defined  limit. 

The  SNR  given  by  Eq.  (78)  corresponds  to  the  SKE/BKE 
decision  task  using  a  discrete-array  detector,  as  Eq.  (13) 
gives  the  SNR  for  the  SKE/BKE  decision  task  for  screen- 
film.  Strictly,  these  formulas  do  not  apply  to  the  task  of 
detection  when  the  observer  does  not  know  the  position  of 
the  object  being  imaged.  For  detecting  a  signal  of  unknown 
location,  one  can  calculate  the  ideal  observer’s  SKE/BKE 
6j{r)  for  each  possible  position  r  of  the  object.  A  common 
strategy  is  then  to  apply  a  threshold  to  0j{r).  Under  the 
assumption  of  Gaussian  statistics  with  complete  knowledge 
except  for  position,  the  likelihood  ratio  computed  by  the 
ideal  observer  uses  0j(r)  *n  a  nonlinear  manner35,39,40  that  is 
sensitive  to  peaks  in  0j(r).  In  either  case,  the  values  of  SNR 
given  by  Eqs.  (13)  and  (78)  are  indicative  of  the  efficacy  of 
the  ideal  observer  in  the  more  general  case  of  the  position 
being  unknown. 

For  the  discrete-array  detector,  however,  the  value  of  the 
SNR  for  the  SKE/BKE  case  will  depend  upon  exactly  where 
the  object  is  relative  to  the  lattice.  While  this  variation  can  be 
significant  (for  example  detectors  could  have  interstitial 
spaces  where  objects  completely  disappear),  the  magnitude 
of  the  effect  decreases  for  objects  large  relative  to  the  lattice 
spacing.  Examples  of  this  for  several  simple  model  detectors 
will  be  given  in  the  next  section.  If  the  variation  in  SNR2 


with  position  is  not  too  great,  then  the  spatially  averaged 
value  of  SNR2  will  be  of  use.1  This  spatial  average  can  be 
computed  exactly  by  noting  that  if  an  object  is  shifted  by  a 
displacement  r,  the  Fourier  transform  is  multiplied  by  e2m{' r 
so  that  in  Eq.  (78)  the  sum  over  elements  of  the  reciprocal 
lattice  becomes 

2 

2  (A/(f+fk))r(f+fk)e27rffk*r  ,  (si) 

fk 

where  a  common  factor  independent  of  k(|e27nf  r|=  1)  has 
been  removed.  In  averaging  over  positions  r  in  Eq.  (78),  the 
denominator  of  the  integrand  does  not  depend  upon  r,  and 
the  numerator  is  the  square  of  the  magnitude  of  a  Fourier 
series  in  r,  so  that  in  integrating  over  r  to  obtain  the  average 
over  all  displacements  one  can  apply  Parseval’s  theorem  to 
obtain 


<SNR2)  =  r2J  j^d2 f 


sfjA/(f+fk)|2|r(f+fk)|2 


W(f) 


=r!|  /  <82) 

where  the  second  step  follows  from  noting  that  the  sum  of 
the  integrals  over  each  unit  cell  is  equivalent  to  the  integral 
over  the  entire  plane. 

As  for  the  OTF  and  NPS,  the  film-screen  result,  Eq.  (13), 
can  be  recovered  from  the  discrete-array  result  [Eq.  (78)]  by 
going  to  the  limit  of  a  sufficiently  fine  lattice,  in  which  case 
the  distance  to  the  first  aliased  frequency  is  so  large  that  only 
the  unaliased  term  contributes  to  Eq.  (78).  Similarly,  for  a 
sufficiently  fine  lattice  all  objects  are  large  relative  to  the 
lattice  spacing,  so  that  SNR2  does  not  vary  appreciably  as  the 
object  is  moved  relative  to  the  lattice  spacing.  These  facts 
prompt  the  identification1  of 

GNEQ(f)  =  r2|r(f)|2d>2/W(f),  (83) 

as  a  generalization  of  the  concept  of  noise  equivalent  quan¬ 
tum  flux  (NEQ),  where  is  the  incident  x-ray  flux,  and 

GDQE(  f)  =  T  2 1  T(  f)  | 2  d>/  W(  f) ,  (84) 

as  a  generalized  detective  quantum  efficiency  (DQE).  These 
results  parallel  the  screen-film  theory,  except  that  factors  of 
fluence  appear  in  the  numerator  as  the  response  of  digital 
detectors  is  linear  with  fluence  [Eq.  (33)]  while  film  density 
is  linear  with  respect  to  the  log  of  fluence  [Eq.  (1)].  While 
Eq.  (82)  is  exact  in  the  context  of  the  assumptions  we  have 
made  about  the  detector,  SNR2  enters  nonlinearly  into  other 
quantities  such  as  the  various  probabilities  of  misclassifica- 
tion  for  a  given  operating  point  (sometimes  called  the  false 
positive  fraction  and  the  false  negative  fraction  in  ROC 
methodology).  However,  when  the  variation  in  SNR2  is  not 
too  large,  perhaps  as  measured  by  the  rms  (root-mean-square 
variation  in  SNR2),  then  the  GNEQ  and  spatially  averaged 
SNR2  can  be  considered  a  useful  summary  of  the  efficacy  of 
the  detector. 

In  this  paper  we  have  applied  the  concept  of  an  ideal 
observer  directly  to  the  digital  data.  The  results  obtained  are 
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implicit  in  the  work  of  Giger  et  aL, 13,25 MA1  but  Giger  et  al. 
concentrates  on  issues  of  display  and  models  of  human  vi¬ 
sual  response  to  the  displayed  data.  As  these  tasks  are  decou¬ 
pled  from  image  acquisition  for  digital  systems,  it  is  worth 
considering  figures  of  merit  for  the  data  acquisition  system 
independent  of  the  display,  as  done  here.  The  results  of  this 
section  also  follow  as  limiting  cases  of  the  work  of  Barrett 
et  al.1  Of  particular  note.  Sec.  VA1  discusses  a  simple  bin¬ 
ning  detector  and  obtains 

"  Ag2 

SNR2=2  — F1.  (85) 

m~  i  <Jm 

where  crm  is  the  uncorrelated  noise  in  the  mth  detector  and 
A  gm  is  the  expected  change  in  the  data  value  at  the  mth 
detector  which  would  be  caused  by  the  signal.  This  particular 
result  can  be  obtained  directly  from  first  principles  based  on 
counting  statistics  in  each  detector  element.  In  the  stationary 
case,  (rm  =  (T  is  a  constant,  so  in  Eq.  (78)  W(f)  =  o-2|A|  and 
the  numerator  [using  Eq.  (38)  and  Parseval’s  identity]  be¬ 
comes  |A|2m|A(D(rm))|2,  again  recovering  the  result  [Eq. 
(85)]  based  on  counting  statistics  for  uncorrelated  noise.  It  is 
worth  noting  that  if  one  does  not  include  the  aliased  terms  in 
the  numerator  of  Eq.  (78)  (perhaps  on  the  grounds  that 
aliased  signals  are  not  useful),  the  value  of  SNR2  will  be 
underestimated.  The  aliased  response  is  part  of  the  physical 
response  of  the  detector,  and  in  this  case  the  aliased  terms 
will  add  coherently  in  such  a  manner  as  to  bring  the  calcu¬ 
lated  value  of  the  SNR2  up  to  the  value  in  Eq.  (85)  obtained 
from  counting  statistics. 


VII.  MODEL  DETECTORS 

To  give  a  feel  for  the  implications  of  the  above  theory,  the 
capabilities  of  detectors  with  reasonably  realistic  parameters 
will  now  be  investigated.  The  modeling  is  somewhat  simplis¬ 
tic,  but  sufficient  to  demonstrate  several  interesting  proper¬ 
ties,  such  as  the  dependence  of  SNR  on  the  position  of  the 
object  being  imaged,  and  certain  trade-offs  inherent  in  such 
detectors,  particularly  those  trade-offs  related  to  the  possible 
suppression  of  input  spatial  frequencies  above  the  frequen¬ 
cies  supported  by  the  lattice.  The  incident  x-ray  fluence 
has  a  white  Wiener  spectrum,  VT/(f)  =  d).  Among  other  sim¬ 
plifications,  which  will  be  discussed  in  more  detail  at  the  end 
of  the  section,  we  assume  100%  of  the  x-rays  interact.  Each 
x-ray  undergoes  a  stochastic  amplification,  characterized  by 
an  average  of  m  secondary  quanta  per  x-ray  with  crm  —  yfm 
for  a  Poisson  process,  and  the  secondary  quanta  undergo  a 
stochastic  scattering  process,  with  a  spread  function  Ps  and 
transfer  function  T$ ,  before  being  “binned”  by  the  detector 
elements.  The  result  is  an  average  of  md>  secondary  quanta 
per  unit  area  on  the  detector  with  a  pre- sampled  Wiener 
spectrum  given  by43 

W,(f)  =  [/n2W,(f)  +  <I>^-/n(&]|ri(f)|2+/n<I>.  (86) 

For  a  square  lattice  with  spacing  L,  binning  can  be  consid¬ 
ered  as  a  deterministic  convolution  with  rect  functions  rep¬ 
resenting  the  detector  regions,  so  that  with 


sin(  irLfx)  \  /  sin(  TrLfy)\ 

ttL/x  J\  ttL/x  /’ 

the  digital  noise  power  spectrum  can  be  written 


(87) 


m©= 4 2  |A|2<D(W2|r,(f+fk)|2+m)|rfe(f+fk)|2 

.  fk 


+  WE,  (88) 

where  the  factor  of  1/m2  is  introduced  so  that  digital  values 
will  correspond  to  x-ray  count  and  WE  is  the  electronic 
noise.  With  the  present  conventions  the  gray-scale  character¬ 
istic  is  set  to  r  =  |  A  | .  A  simplification  can  be  achieved44 
using 


sin(7r(jt  +  w))\ 2 
7 r(x  +  /2)  / 


(89) 


for  any  x,  which  can  be  proven  by  applying  Parse val's  theo¬ 
rem  to  the  Fourier  series  for  e2mxy  for  y  e  [  —  0.5, 0.5].  The 
experimentally  observable  transfer  function  (as  obtained,  for 
example,  by  the  slanted-edge  technique,  cf.  Sec.  IV)  contains 
the  effects  of  stochastic  scatter  and  binning,  thus  T{ f) 
=  Ts(f)Tb( f),  so  that 

W(f)/|A|  =  d)|A|2  |r(f+fk)|2+— +  w£/|a|  (90) 

fk  m 

is  the  Wiener  spectrum  of  the  model  detector,  with  the  aver¬ 
age  number  of  x-rays  per  pixel  being  <F|A|.  The  summation 
over  aliases  in  Eq.  (88)  is  often  referred  to  as  “noise  alias¬ 
ing.”  The  division  into  aliased  and  unaliased  components  is 
useful  for  modeling  a  variety  of  detectors,  but  it  should  be 
noted  that  this  division  is  generally  not  directly  experimen¬ 
tally  accessible,  at  least  not  without  modifying  the  detectors, 
and  that  in  principle  there  could  be  devices  which  are  sta¬ 
tionary,  and  therefore  have  Wiener  spectra,  but  for  which  the 
division  of  the  NPS  into  aliased  and  unaliased  components  is 
not  useful. 

It  is  useful  to  choose  values  of  the  parameters  in  the 
model  which  are  representative  of  detectors  of  current  clini¬ 
cal  interest,  as  this  can  help  in  the  understanding  of  the  phys¬ 
ics  which  determines  the  performance  of  these  devices,  but 
detailed  modeling  for  quantitative  comparison  to  actual  de¬ 
vices  is  beyond  the  scope  of  this  article.  We  assume  a  square 
lattice  with  spacing  of  L~ 0.143  mm,  operation  at  an  expo¬ 
sure  corresponding  to  |A|<F=  1400  x-rays  per  pixel,  and  an 
amplification  factor  of  m  =  1000.  For  the  stochastic  transfer 
function  Ts  we  consider  three  possibilities:  a  “blur-free” 
detector  for  which  Ts(f)~  1,  typical  of  photoconductive 
arrays,45  and  two  “alias-free”  detectors  whose  stochastic 
transfer  functions  are  of  the  form  T5(f)  =  ^“x^  with 
\ =0.463  and  \=0.34  mm,  which  approximates  the  transfer 
function  for  evaporated  Csl.46-48  Typically  electronic  noise 
V(  W£)/|A|  is  on  the  order  of  3-5  x-rays,  so  values  of  0,  42, 
and  82  cover  the  range  of  values  for  $>\A\/m  +  WE/\A\. 

The  transfer  functions  for  these  models  are  shown  in 
Fig.  2.  As  the  pixels  are  symmetric  with  respect  to  inversion 
through  their  centers  [i.e.,  for  the  PSF,  P(r)  =  P(-r),  and 
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<D|A|/m  +  We/|A|  =  0 


Fig.  2.  The  optical  transfer  functions  of  three  model  detectors.  The  “blur- 
free"  detector  bins  the  secondary  quanta  without  smoothing,  while  for  the 
“alias-free"  detectors  the  distribution  of  secondary  quanta  is  smoothed  by 
an  exponential  MTF  ( e~K |f|)  before  binning.  Data  are  shown  as  a  function  of 
the  magnitude  of  the  spatial  frequency  for  several  angles. 


P(r)  is  a  real  number],  the  imaginary  part  of  the  transfer 
function  is  identically  zero,  so  only  the  real  part  need  be 
graphed.  The  OTF  is,  of  course,  a  function  of  two  variables, 
fx  and  fy .  To  show  this,  we  plot  the  OTF  as  a  function  of  the 
magnitude  of  the  frequency  vector  for  three  angles  relative  to 
an  axis  of  the  detector.  For  the  blur-free  detector,  the  transfer 
function  is  simply  the  product  of  the  sines  in  the  two  direc¬ 
tions  induced  by  the  binning  operation.  The  OTF  of  the  blur- 
free  detector  is  nonzero  well  beyond  the  highest  frequency 
supported  by  the  lattice.  Any  component  of  an  input  signal  at 
these  higher  frequencies  will  contribute  to  a  lower  frequency 
alias  in  the  output,  as  per  Eq.  (38),  and  while  it  is  not  obvi¬ 
ous  from  the  point  of  view  of  frequency  space  the  sum  over 
aliases  in  Eq.  (38)  will  be  precisely  equivalent  to  the  detector 


Wiener  Spectra 


Fig.  3.  The  Wiener  spectra  W(f)/|A|  for  the  three  model  detectors  as  in  Fig. 
2.  The  residual  additive  white  noise  4>|A|/m  +  lV£/|A|  has  been  set  to  0. 


Fig.  4.  GDQE  as  a  function  of  frequency  for  the  three  model  detectors  as  in 
Fig.  2,  with  the  residual  white  noise  4>|A|/m  +  We/|A|  set  to  0. 


simply  binning  each  incident  x-ray.  For  the  alias-free  detec¬ 
tors,  there  is  very  little  response  to  frequencies  beyond  those 
supported  by  the  lattice.  For  each  detector,  three  angles  are 
plotted,  but  the  angular  dependence  for  the  alias-free  detec¬ 
tors  is  small  enough  to  not  be  apparent  on  the  graph. 

The  Wiener  spectra  are  shown  in  Fig.  3.  Again,  instead  of 
plotting  a  function  of  two  variables,  f x  and  /v ,  we  plot  the 
NPS  as  a  function  of  the  magnitude  of  the  frequency  vector 
for  three  angles,  0=0,  27,  and  45°  (27°  corresponds  to  a 
slope  of  1:2  relative  to  the  lattice).  However,  the  NPS  shows 
little  angular  dependence.  For  the  Wiener  spectrum  one  only 
needs  to  look  at  frequency  values  supported  by  the  lattice, 
i.e.,  fxe[-  1/2L,1/2L]  and  fy  e [ -  1/2L,1/2L]  (for  conve¬ 
nience  one  can  consider  W  to  be  periodic  in  the  frequency 
plane).  Thus,  at  0=0°  one  only  needs  to  graph  up  to  1/2 L 
=  3.5  mm-1,  but  at  0=45°  the  frequencies  are  on  the  diago¬ 
nal  of  the  square,  so  one  goes  up  to  \j2/2L — 4.9  mm  1 .  At 
0=27°,  one  goes  up  to  l/(2Lcos  0)=4mm_1.  For  this  graph 
the  constant  offset  <I>|A|/m  +  W£-/|A|  has  been  set  to  zero. 
For  the  blur-free  detector,  the  NPS  is  flat,  which  follows 
mathematically  from  Eq.  (89)  and  the  fact  that  F(f)  for  these 
detectors  is  simply  related  to  sine  functions,  or  more  physi¬ 
cally  by  noting  that  for  a  detector  which  simply  bins  incident 
x-rays  adjacent  cells  will  be  uncorrelated  so  the  NPS  is  flat. 
For  the  alias-free  detectors,  the  NPS  is  suppressed  by  factors 
of  the  square  of  the  transfer  function. 

The  GDQE  as  a  function  of  frequency  are  shown  in 
Figs.  4-6  for  a  range  of  values  of  the  residual  white  noise 
<b\A\lm  +  WEl\A\.  For  each  graph,  values  are  plotted  as  a 
function  of  the  magnitude  of  the  spatial  frequency  for  angles 
0,  27,  and  45°  relative  to  an  axis.  In  each  case,  the  GDQE 
falls  off  most  quickly  at  0=0°  and  least  quickly  at  0=45°, 
which  represents  the  fact  that  on  the  diagonal  the  sampling 
rate  is  increased  by  a  factor  of  V2.  In  the  case  where  the 
residual  white  noise  is  zero  (Fig.  4),  the  GDQE  of  the  blur- 
free  detector  drops  like  the  square  of  a  sine  function.  For  the 
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SNR2  as  a  function  of  displacement 


{50p.  wire  20  mm  long  parallel  io  axis) 


Fig.  5.  GDQE  as  a  function  of  frequency  for  the  three  model  detectors  as  in 
Fig.  2,  with  the  residual  white  noise  <$>\A\lm  +  WEl\A\  set  to  42. 


alias-free  detectors,  the  GDQE  remains  at  nearly  unity  up 
to  the  lattice  cutoff,  the  factor  of  T2(  f)  canceling  the  same 
factor  in  the  colored  part  of  the  noise.  The  GDQE  of  the 
blur-free  detector  shows  some  response  beyond  the  lattice 
cutoff.  Though  small,  this  portion  of  the  GDQE  is  physical 
and  it  will  be  shown  that  the  responses  to  the  aliased 
frequencies  can  not  be  trivially  dismissed.  With  the  addition 
of  residual  white  noise  the  GDQE  of  all  three  models  is 
reduced,  as  shown  in  Fig.  5  (<t>|A|/m  +  We/\A  |  =  42)  and 
Fig.  6  (0>\A\lm  +  WEl\A\  =  %2),  These  figures  illustrate  that 
the  alias-free  detectors  are  more  sensitive  to  sources  of  re¬ 
sidual  white  noise  than  blur-free  detectors.  Indeed,  in  Fig.  6 
the  blur-free  detector  now  has  higher  GDQE  than  the 
\ =0.463  mm  detector  even  at  low  spatial  frequencies.  From 
the  spatial  point  of  view  this  is  quite  reasonable.  A  detector 


GDQE 


<t>|A|/m  +  W£/|A|  =  82 


Fig.  6.  GDQE  as  a  function  of  frequency  for  the  three  model  detectors  as  in 
Fig.  2.  with  the  residual  white  noise  4>|A|/w  + W£/|A|  set  to  82. 


Fig.  7.  Relative  SNR2  as  a  function  of  position  for  a  50  jtim  wide,  20  mm 
long  wirelike  object  parallel  to  one  axis  of  the  detector.  The  horizontal  axis 
of  the  graph  gives  the  displacement  of  the  wire,  so  that  at  0  mm  the  wire  is 
over  a  single  column  of  sensitive  elements,  while  at  0.07  mm  the  wire 
straddles  two  columns.  The  verticle  scale  is  arbitrary  (dependent  upon  the 
contrast  of  the  wire). 


whose  transfer  function  is  designed  to  remove  aliases  has  a 
relatively  wide  point  spread  function,  and  therefore  a  rela¬ 
tively  wide  autocovariance  function.  The  ideal  observer 
makes  use  of  digital  values  in  array  elements  whose  distance 
from  the  position  of  the  signal  is  up  to  several  times  the 
lengths  of  the  tails  of  these  functions,  so  for  the  same  input 
signal  the  ideal  observer  will  have  to  integrate  over  a  larger 
region  on  an  alias-free  detector  and  thus  will  be  more  sensi¬ 
tive  to  any  residual  uncolored  noise. 

While  GDQE  is  directly  related  to  the  average  value  of 
SNR2  by  Eq.  (82),  high-frequency  signals  (i.e.,  x-ray  pro¬ 
jection  images  of  small  objects  or  objects  whose  projected 
density  varies  quickly  with  position)  can  demonstrate  signifi¬ 
cant  changes  in  SNR2  with  position.  To  explore  this,  con¬ 
sider  the  SKE/BKE  task  associated  with  an  object  50  /im 
wide  and  20  mm  long.  Wires  of  this  width  have  been  used  in 
neurological  and  cardiovascular  stents.49  Figure  7  shows  the 
SNR2  for  such  an  object,  oriented  parallel  to  an  axis  of  the 
detector,  as  a  function  of  displacement  in  the  direction  of  the 
shorter  (50  /mi)  axis.  The  scale  of  the  vertical  axis  is  arbi¬ 
trary  as  we  won't  set  the  inherent  contrast  of  the  signal.  In 
Fig.  7,  the  0  mm  displacement  corresponds  to  the  signal 
being  centered  over  the  sensitive  elements  of  the  detectors. 
For  the  blur-free  detector,  at  0  mm  displacement  the  signal 
falls  into  a  single  column  of  detectors,  so  the  SNR2  corre¬ 
sponds  to  counting  statistics  for  one  column  of  elements  20 
mm  long.  The  SNR2  is  constant  until  the  0.05  mm  mark, 
after  which  the  signal  is  shared  between  two  columns  of 
detectors,  resulting  in  a  drop  in  SNR2.  Physically,  the  num¬ 
ber  of  x-rays  attenuated  by  the  object  is  independent  of  its 
position,  but  at  a  displacement  of  0.7  mm  the  signal  is  shared 
equally  between  two  columns,  and  since  the  noise  is  as¬ 
sumed  to  be  uncorrelated,  the  variance  in  the  total  counts  in 
the  two  columns  is  twice  that  of  one  column,  so  the  SNR2  is 


Medical  Physics,  Vol.  27,  No.  10,  October  2000 


2432 


M.  Albert  and  A.  Maidment:  Linear  response  theory 


Table  I.  SNR2  averaged  over  position  and  orientation  for  the  projection 
image  of  an  object  0.05  by  20  mm.  As  the  SNR2  scales  with  the  square 
of  the  contrast  of  the  image,  only  relative  values  are  meaningful.  The  ± 
represent  the  rms  (root-mean-square)  fluctuations  in  the  SNR2  with  position, 
not  the  statistical  uncertainties. 


SNR2  for  50  fim  wide,  20  mm  long  object 

<D|A|/m  +  W£/|A|=0 

=42 

=  82 

Blur-free 

1.81  ±0.02 

1.79±0.02 

1.73  ±0.02 

Alias-free  (\= 0.463  mm) 

2.16±0.16 

1.87  ±0.07 

1.48  ±0.02 

Alias-free  (\=0.34  mm) 

2.15±0.15 

1.99  ±0.10 

1.68  ±0.05 

reduced  by  half.  The  alias-free  detectors  show  less  sensitivity 
to  position,  as  the  signal  is  always  shared  between  multiple 
columns.  As  before,  curves  are  shown  for  three  values  of  the 
residual  white  noise  <&|A|/m  +  WEl\A\(02A2y  and  8“).  For 
all  models,  the  SNR2  drops  as  the  residual  white  noise  in¬ 
creases,  but  this  effect  is  greater  for  the  alias-free  models. 
Table  I  gives  the  average  SNR2  for  detection  of  the  0.05 
wide,  20  mm  long  wire,  now  averaged  over  both  position 
and  orientation.  Additionally,  the  root-mean-square  variation 
in  SNR2  is  given,  to  indicate  the  degree  to  which  the  detect¬ 
ability  of  the  wire  would  vary.  Again,  the  alias-free  detectors 
give  a  higher  SNR2  if  the  residual  white  noise  is  kept  suffi¬ 
ciently  low.  In  calculating  the  SNR2  of  the  projection  of  the 
wire  using  Eq.  (78),  if  the  summation  over  aliases  is 
dropped,  the  resulting  integral  decreases  by  about  10%.  Thus 
the  contributions  of  the  aliased  signal  to  the  SNR  are  not 
always  negligible. 

Somewhat  speculatively  one  can  consider  tasks  which  de¬ 
pend  upon  higher  frequency  components  of  the  signal. 
Consider  a  5  mm  square  with  10%  contrast,  and  a  second 
square  whose  edges  have  been  smoothed  by  convolving  with 
a  0.15  mm  rect  function,  so  that  the  resulting  signal  “ramps 
up'  ’  over  a  distance  of  0.3  mm.  The  SNR2  for  the  SKE/BKE 
task  of  distinguishing  between  these  two  objects  is  given  in 
Table  II.  It  is  interesting  to  note  that,  mathematically,  the 
SNR2  is  sufficiently  large  that  the  ideal  observer  can  perform 
this  task  efficiently,  although  whether  a  human  could  do  this 
is  questionable.  On  the  other  hand,  edge  detection  is  impor¬ 
tant  both  computationally  and  probably  as  part  of  the  strat¬ 
egy  of  human  observers,  so  the  ability  to  perform  this  task  is 
not  a  priori  irrelevant.  Again,  the  alias-free  detectors  have  a 
higher  SNR2  if  the  residual  white  noise  is  zero,  but  as  the 
task  now  depends  more  heavily  on  the  higher  frequency 

Table  II.  SNR2  averaged  over  position  and  orientation  for  distinguishing 
between  the  projection  of  a  5  mm  square  with  sharp  boundaries  and  a  square 
whose  boundaries  ramp  up  over  a  0.3  mm  region.  Normalization  corre¬ 
sponds  to  a  10%  contrast  with  an  x-ray  flux  corresponding  to  1400  x-rays/ 
pixel.  The  ±  represents  the  rms  fluctuations  in  the  SNR^  with  position,  not 
the  statistical  uncertainties. 


SNR2  for  discontinuity  vs  slope 


<D|A|/m  +  W£/|A|=0 

=  42 

=  82 

Blur-free 

14±2 

14±2 

13±2 

Alias-free  (X =0.463  mm) 

16±3 

12±1 

7.6±0.5 

Alias-free  (X=0.34  mm) 

16±3 

14±2 

10±1 
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Table  III.  SNR2  averaged  over  position  and  orientation  for  detecting  the 
projection  of  a  0.5  mm  square.  Normalization  is  arbitrary.  The  ±  represents 
the  rms  fluctuations  in  the  SNR2  with  position,  not  the  statistical  uncertain¬ 
ties. 


SNR2  for  0.5  mm  square 

Blur-free 

Alias-free 

<t>|A|/m  +  WE/|A|  =  0 
1.91  ±0.01 

2.08  ±0.01 

=42 

1.89±0.01 

1.99±0.01 

=  82 

1.83±0.01 

1.789±0.003 

(X =0.463  mm) 
Alias-free 

2.07±0.01 

2.01  ±0.01 

1.855  +0.006 

(X=0.34  mm) 

components  of  the  signal,  the  alias-free  detectors  are  more 
sensitive  to  residual  white  noise,  with  the  crossover  at 
<F|A|/m  + W£/|A|==42.  The  detection  of  a  0.5  mm  square  is 
shown  in  Table  III.  Here,  the  lower  (but  nonzero)  frequen¬ 
cies  dominate  the  response  of  the  detector,  so  that  in  general 
the  antialiasing  detectors  gain  from  the  removal  of  the 
aliased  noise  without  losing  any  signal. 

It  is  interesting  to  note  that  the  SNRs  for  the  tasks  and 
models  described  above  do  not  vary  greatly.  Many  factors 
not  considered  here  will  greatly  effect  the  performance  of 
real  detectors,  beginning  with  the  fact  that  less  than  100%  of 
the  incident  x-rays  will  produce  secondary  quanta.  The  color 
of  the  Wiener  spectrum  need  not  be  the  same  as  the  transfer 
function,  due  to,  for  example,  x-rays  interacting  at  various 
depths  in  the  detector.51  The  efficiency  of  collection  of  the 
secondary  quanta  can  also  have  significant  effects.52  For  Csl 
detectors,  of  which  our  “alias-free”  detector  is  a  rough 
model,  the  fill  factor  is  a  minor  effect  so  long  as  the  ampli¬ 
fication  m  is  sufficiently  large.48  For  selenium  detectors,  of 
which  our  “blur-free”  detector  is  an  approximation,  it  is 
possible  to  have  an  effective  fill  factor  significantly  greater 
than  the  geometric  fill  factor  of  the  TFT  array.53,54  In  any 
case,  our  purpose  here  is  merely  to  indicate  some  of  the 
issues  which  must  be  faced  in  quantifying  digital  detectors  of 
these  types.  In  addition,  we  did  not  consider  geometric  fac¬ 
tors  such  as  x-ray  focal  spot  size  and  x-ray  parallax  which 
reduce  the  high-frequency  content  of  the  incoming  signals. 

VIII.  DISCUSSION  AND  CONCLUSION 

The  results  of  this  paper  set  up  a  framework  for  quantita¬ 
tive  measurements  of  digital  systems  in  a  manner  analogous 
to  the  now  common  analysis  of  screen-film  systems  in  terms 
of  the  gray-scale  transfer  characteristic,  the  optical  transfer 
function,  the  Wiener  spectrum,  and  signal-to-noise  ratio.  The 
logic  of  this  framework  is  by  design  close  to  the  classic  work 
on  film-screen  systems.  While  many  pieces  of  this  argument 
have  appeared  in  the  work  of  others,  as  noted  throughout  the 
text,  it  seemed  desirable  to  produce  a  coherent  systematic 
exposition.  The  results  can  be  seen  as  appropriate  limits  of 
those  of  Barrett  et  al ,  but  the  theoretical  construction  here 
emphasizes  the  parallels  with  the  classic  results  on  screen- 
film  systems.  While  detectors  consisting  of  discrete  elements 
do  not  have  continuous  translational  symmetry,  the  remain¬ 
ing  discrete  symmetry  allows  one  to  use  the  appropriate  Fou¬ 
rier  technique.  The  advantage  of  this,  as  in  the  screen-film 
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case,  is  that  one  can  explicitly  solve  for  the  mask  function  of 
the  ideal  observer  [see  Eq.  (75)]  and  thereby  obtain  the  op¬ 
timal  SNR.  As  for  screen-film,  this  formula  can  be  inter¬ 
preted  as  the  ratio  of  the  square  of  the  signal  in  each  fre¬ 
quency  bin  to  the  noise  in  each  frequency  bin,  as  measured 
by  the  Wiener  spectrum,  integrated  over  bins. 

We  have  investigated  several  highly  idealized,  but  not 
completely  unrealistic,  models  of  detectors,  and  illustrated 
some  of  the  issues  inherent  in  various  design  decisions.  This 
analysis  is  incomplete  and  intended  to  point  toward  issues 
which  could  be  addressed  in  other  work.  However,  our  re¬ 
sults  suggest  that  for  typical  tasks  the  detectability  of  objects 
as  determined  by  SNR2  is  not  drastically  affected  by  the 
decision,  in  and  of  itself,  to  suppress  or  not  to  suppress 
aliases.  In  any  real  device,  of  course,  this  design  decision  is 
linked  to  many  other  parameters.  This  article  should  be  of 
use  in  clarifying  what  is  actually  experimentally  measured  in 
testing  such  devices. 

The  results  presented  here  are  exact  for  the  SKE/BKE 
task  as  approached  by  the  ideal  observer  under  the  assump¬ 
tions  of  linearity,  homogeneity,  and  stationarity.  However, 
each  of  these  assumptions  is  only  approximately  true  in  prac¬ 
tice.  The  finite  extent  of  real  detectors  trivially  shows  that 
they  are  not  homogeneous,  but  for  a  variety  of  tasks  edge 
effects  are  negligible.  More  importantly,  many  digital  detec¬ 
tors  in  practice  show  significant  inhomogeneity  and  nonsta- 
tionarity.  The  work  of  Barrett  et  ai  is  sufficiently  general  to 
cover  these  cases.  Further,  the  inhomogeneity  and  nonsta- 
tionarity  of  a  given  instrument  often  occur  in  ways  which  are 
different  for  each  individual  device,  so  that  while  the  extra 
information  is  relevant  to  the  particular  device  one  has  mea¬ 
sured,  the  extra  information  is  often  not  generalizable  to 
other  devices  of  the  same  manufacture.  This  extra  informa¬ 
tion  is  useful  for  optimizing  certain  tasks  using  the  particular 
device,  but  of  less  use  in  understanding  a  class  of  devices. 
Additionally,  while  the  signal  detection  task  of  the  ideal  ob¬ 
server  under  SKE/BKE  conditions  certainly  shares  some  fea¬ 
tures  with  the  task  which  human  observers  face,  and  has 
under  many  conditions  been  shown  to  correlate  well  with  the 
ability  of  human  observers  (for  example,  screen-film  images 
of  nylon  beads55),  it  is  still  a  very  idealized  task.  For  ex¬ 
ample,  if  edge-detection  is  important,  then  higher  frequency 
parts  of  the  incoming  signal  become  more  important  than 
would  be  expected  given  simply  the  signal  detection  task. 
From  the  point  of  view  of  a  radiologist,  a  clear  edge  might  be 
used  to  identify  and  distinguish  the  existence  of  a  lesion 
from  a  variation  in  projected  density  of  the  underlying  organ, 
particularly  in  the  presence  of  the  “structured  noise”  of 
other  anatomical  features. 

Clearly  these  are  issues  for  further  study,  but  while  OTF, 
W,  SNR,  and  GDQE  are  certainly  useful  because  they  are 
objectively  measurable  and  in  a  mathematically  precise  man¬ 
ner  are  related  to  tasks  which  approximate  those  of  the  hu¬ 
man  observer,  it  is  worth  remembering  that  measurement  of 
these  quantities  does  not  obviate  the  need  for  observer  stud¬ 
ies,  particularly  with  practicing  radiologists. 
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